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Abstract 

Domain Decomposition Methods for the Electromagnetic Simulation of Scattering from 
Three-Dimensional Structures with Applications in Lithography 

by 

Konstantinos Adam 

Doctor of Philosophy in Electrical Engineering 

University of California, Berkeley 

Professor Andrew R. Neureuther, Chair 

An integrated methodology has been developed for the computer simulation and 
modeling of the electromagnetic scattering from large, non-periodic, two-dimensional lay- 
outs of advanced photomasks (masks with optical proximity correction and phase shifting 
masks). The name domain decomposition method (DDM) was used, since it describes the 
central mechanism of the method. 

Domain decomposition consists of three important steps: First, by virtue of the lin- 
earity of the Kirchhoff-Fresnel diffraction integral, the mask layout is decomposed into a 
set of constituent single-opening masks. Secondly, the rigorous electromagnetic simulation 
of each three-dimensional structure from the set of these single-opening masks is circum- 
vented, and instead, the result for the scattered field is synthesized based on two two- 
dimensional rigorous electromagnetic simulations that model the mask geometry in two 
cross-sectional planes. Subsequently, based on the results of the electromagnetic scattering 
from these two-dimensional geometries, compact equivalent source models are used to 
describe the scattered fields on a reference plane. These models are constructed in such a 
way as to minimize the error in the part of the diffraction spectrum that is passing through 
the projection system allowing accurate and efficient image simulation. 

Excellent accuracy in the calculation of the near scattered fields of better than 99% 
(in a normalized mean square error sense) compared with the fully rigorous mask model 
has been achieved, accompanied by speed-up factors for the total simulation time in excess 
of 200. A further revision of the method consisting of the decomposition of the layout into 
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edges (edge-DDM) allows for easier algorithmic implementation. The algorithm of the 
edge-DDM was programed in the MATLAB environment and together with TEMPEST 
cross-sectional simulations resulted in speed-up factors for the total simulation time of 
172,800 (lsec. vs. 2days) using a library of pre-calculated edge-diffraction simulations. 
The normalized mean square error of the near field results between the edge-DDM and full 
three-dimensional simulation is less than 1%. A 12um by 16nm layout of a three level alter- 
nating PSM, that is out of reach for fully rigorous methods, was simulated in under lmin. 

The domain decomposition method was extended in two ways: Energy cross-cou- 
pling between neighboring apertures in alternating phase shift masks was modeled through 
simulation, enabling the accurate modeling of masks with large vertical topography. 
Another revision of the domain decomposition method suitable for handling the case of 
scattering from masks when phase defects are present was also realized, enabling rapid 
defect printability assessment. 

Finally, useful engineering design data relevant to the design of optical proximity 
correction were produced and the electromagnetic behavior of isolated phase defects was 
examined to understand their interaction with layout features. 



Professor A.R. Neureuther 
Committee Chairman 
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^ Introduction 



ONE of the most extraordinary achievements of engineering ingenuity has been the 
ability of the semiconductor industry to make a law out of a prediction in 1965 by Gordon 
Moore, meant to be valid for no more than 10 years, and according to which "the complex- 
ity of an integrated circuit, measured by the number of components (transistors, capacitors, 
diodes, resistors) per chip, would double approximately every two years" [58], More than 
thirty five years later and after an increase in complexity approaching a factor of one hun- 
dred million (the one billion transistor integrated circuit is expected to appear before the 
end of the decade), the end of this spectacular exponential growth has yet to happen. Argu- 
ably, by far the biggest contributor to this progress has been what is known as optical pro- 
jection lithography. This term loosely implies the physical process by which the integrated 
circuit is printed on a semiconductor wafer while all the circuitry information travels 
through a complex optical system in the form of ultraviolet light rays. 

The end of optical lithography has been unsuccessfully predicted many times over 
the past ten years. However, everybody now agrees that certain physical and excessive cost 
barriers will eventually, probably in this decade, necessitate the replacement of this tech- 
nology. Extreme ultraviolet lithography (EUVL) has emerged as the most prominent 

replacement candidate 1 . Although EUV lithography is a profoundly different technology, 
the main physical principles are nevertheless the same, since it will again be through infor- 
mation carried by higher frequency photons that the printing of integrated circuits will be 
realized. 



1. See for example: "Extreme Ultraviolet Lithography: Will it be Ready in Time?" in November 2001 issue 
of IEEE Spectrum. 



1 



c 



c 



Optical photolithography will not go out without a fight. This fight is currently 
taking place in various fronts, most of which are grouped under what is known as resolution 
enhancement techniques (RET). RET is nothing more than an elaborate way to convey the 
"tricks" that engineers are using to make further increases in integrated circuit complexity 
possible. Two important such tricks are optical proximity correction (OPC) and phase shift 
masks (PSM), that are central to this work. 

The widely accepted practice of using physical modeling through computer simu- 
lation in order to predict, design and model the outcomes of the lithography process was 
used in this research. Such a practice is justified for various reasons having to do with the 
reduced cost of conducting "virtual experiments" through computer simulations as opposed 
to costly real experiments in state-of-the-art integrated circuit fabrication facilities, and also 
with the fact that the various physical phenomena can be easily decoupled and better under- 
stood. However, the value of modeling through computer simulation goes far beyond that. 
Certain resolution enhancement techniques, such as the optical proximity correction, the 
optimization of phase-well topography in alternating phase shift masks and the concurrent 
illumination/mask optimization rely heavily on modeling through computer simulation 
[24]-[27], [63]-[64], [76]. 

This research focuses on the modeling of optical phenomena (scattering and prop- 
agation) that relate to photomasks. In particular, the techniques developed aim to bridge the 
gap between fully rigorous modeling methods that are not applicable in larger scale prob- 
lems because of enormous computational requirements (speed and memory) and faster but 
less accurate methods that derive from light diffraction theory. 

1.1. Organization of text 

The operation of an optical projection system is introduced in Chapter 2 with spe- 
cial focus on the details that are relevant to this work. Elements of the scalar diffraction 
theory for partially coherent light are reviewed and the application of rigorous models for 
the mask diffraction problem are more extensively examined. In particular, issues that 
relate to the modeling of an extended illumination source and the origin of polarization 
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effects and their simulation are discussed in detail. Various simulation approaches for the 
mask diffraction problem are briefly reviewed and categorized. 

Chapter 3 deals with optical proximity correction from a design point-of-view. The 
purpose is to uncover the type and magnitude of the corrections that are influenced by the 
geometrical parameters of assist features, such as scatter bars and corner serifs. The need 
for rigorous mask diffraction simulation owing to the small (wavelength-sized) size of 
assist features is critically examined and a simple perturbation approach is used to adapt 
the scalar diffraction model to the correct solution in the case of binary Chromium-on-glass 
(COG) masks. 

The body of the more significant work begins with Chapter 4. A methodology is 
systematically developed that permits accurate electromagnetic simulation of the problem 
of mask diffraction for masks with large and arbitrary two-dimensional layouts. The start- 
ing point is the linearity of the Kirchhoff-Fresnel diffraction integral and what it implies 
for complex two-dimensional diffracting screens. Based on this linearity, a decomposition 
of the domain into a set of constituent parts is performed that is in principle the same as the 
field stitching method of Layet et al [49] and Prather et al. [72] for solving the problem of 
large, one-dimensional diffraction gratings. Following this first decomposition, a second 
division of the elemental two-dimensional problems into a set of one-dimensional diffrac- 
tion problems takes place. The accuracy of this step is explored and the physical reasons 
behind this accuracy are pointed out. Specifically, corner effects and polarizations perpen- 
dicular to the excited one are shown to be insignificant. Next, the act of the optical system 
as a low pass filter on the diffraction spectrum is pointed out and the idea of spectrally 
matched simpler models is developed. The algorithmic implementation of the method is 
also outlined. 

In Chapter 5 the problem of energy cross-coupling between neighboring apertures 
of an alternating phase shift mask is treated. This is a critical step for the practical imple- 
mentation of the domain decomposition framework, since the cross-communication, if not 
modeled, can severely confound the accuracy. From a single scattering simulation the 
physical mechanism of cross-talk becomes clear and a model is built to quantitatively take 
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cross-talk into account. Details that pertain to the dependence of cross-talk on the geomet- 
rical characteristics of the mask, the insertion point when the model is necessary and ideas 
for cross-talk elimination are also given. 

The decomposition methods of Chapters 4 and 5 are further expanded in Chapter 6. 
Here, one extra division of the elemental one-dimensional diffraction problems into a set 
of edges takes place. The idea is similar to the technique of using edge-diffraction informa- 
tion for the field reconstruction of one-dimensional dielectric gratings by Kettunen et al. 
[48], but important differences are pointed out. The critical step of edge decomposition 
brings tremendous versatility and the algorithmic implementation is now less complicated. 
The accuracy of the method is, again, systematically examined. The implementation of the 
proposed algorithm into the MATLAB environment and examples of large, arbitrary lay- 
outs that are out of reach for rigorous methods are presented. The idea of spectrally matched 
elements that are compact is also used for the edge diffraction problem. 

Phase defects in alternating phase shifting masks are the topic of Chapter 7. In the 
first part, the electromagnetic fingerprint of phase defects based on their type (glass protru- 
sions or cavities) and size is quantified. Fundamental differences between the electromag- 
netic scattering of the two types of defects are revealed and used to explain the differences 
of linewidth variations when identically sized defects are present. Recent publications have 
provided experimental support of the results that were first obtained through computer sim- 
ulation [93], [23]. In the second part of the Chapter, the domain decomposition method is 
adapted for simulation of defective alternating phase shifting masks in a way that should be 
suitable either for rapid defect printability assessment or die-to-database comparisons in 
mask inspection systems. 

Finally, a summary of the most important results and ideas for future research on 
the methods developed can be found in Chapter 8. 

The organization of each Chapter (other than Chapter 2, which can be considered 
introductory material in its entirety) is such that a short introduction relevant to the subse- 
quent topics is given first, followed by the body of main ideas and concluded with impor- 
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tant results and observations. References to each topic can be found in the beginning of 
each Chapter, but are also scattered throughout, whenever it is appropriate. 

1.2. Thesis contributions 

Chapters 3 and 7 contain a rich amount of engineering data relevant to the design 
of optical proximity correction and phase defects respectively. The successful application 
of a perturbation model to adapt the accuracy of the scalar method is probably the most sig- 
nificant contribution of Chapter 3, whereas the quantitative characterization of phase 
defects and, in particular, the observation that the actual phase shift that they induce 
locally, depends strongly on the lateral dimension of the defect is the highlight of the first 
part of Chapter 7. 

By far though, the most important contribution of this thesis has to be the develop- 
ment of the simulation methodology in Chapters 4-7 based on decomposition of the orig- 
inal mask geometry, that enables rapid diffraction analysis of large masks with arbitrary 
layouts, without sacrificing accuracy. This methodology should have a direct impact in 
model-based optical proximity correction and in the inspection and printability assessment 
of phase defects, where speed is critical. 
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V Tn * photolithography process and its 
simulation 



THE purpose of this Chapter is to provide the theoreticaJ and other background on 
the subsequent thesis material. An effort is made to keep the discussion short, yet concise, 
but some important topics are presented with greater detail. 

2.1. Operation principle of optical photolithography 

The general components of an optical lithography tool, shown schematically in the 
diagram of Figure 2-1, are the illumination system, the projection system, the photomask 
( also called reticle), and the photoresist spun on top of a semiconductor wafer. The opera- 
tion principle of the system is based on the ability of the resist to record an image of the 
pattern to be printed. The mask, already carrying this pattern, is flooded with light and the 
projector forms an image of all mask patterns simultaneously onto (and into) the resist. The 
inherent parallelism of this process is the main reason why optical photolithography is 
favored over any other lithography, since it facilitates a very high throughput of 30-120 
wafers per hour. The light intensity distribution on top of the resist surface is commonly 
referred to as aerial image. The resist itself is a photosensitive material whose chemical 
composition changes during light exposure. The pattern is thereby stored in form of a latent 
(bulk) image within the resist. After exposure has occurred, the resist is developed by 
means of a chemical process that resembles the process of developing photographic film. 
After development, the exposed parts of the resist remain or dissolve depending on its 
polarity (negative or positive respectively). The end-result of the lithography process is a 
more or less exact (scaled or not) replica of the mask pattern on the wafer surface that will 
play the role of a local protective layer (mask) for subsequent processing steps (etching, 
deposition, implantation). 
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Figure 2-1. General diagram of an optical photolithography tool 



Brief descriptions of the illumination system, the projection system and the photo- 
mask are given next. The complex role of the photoresist in the photolithography process 
is beyond the scope of this thesis, and the details of the resist chemistry will not be dis- 
cussed in the following. 

2.1.1. Illumination system 

The role of the illumination system is to deliver a light beam that uniformly trans- 
illuminates the entire reticle. It typically consists of various optical elements, such as 
lenses, apertures, filters and mirrors. The light source is responsible for generating very 
powerful and monochromatic radiation. Power is necessary because it is directly related to 
throughput. Monochromaticity is important because high quality refractive (or reflective 
in the case of EUV lithography) optics can only be fabricated for a very narrow illumina- 
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tion bandwidth. State-of-the-art optical lithography tools employ excimer 1 lasers as their 
light source. Deep Ultra Violet (DUV) lithography is the term used for lithography systems 
with illumination wavelengths A=248nm (excimer laser with KrF), X=193nm (excimer 
laser with ArF) and X=157nm (excimer laser with F 2 ). The successful development of cur- 
rent and future optical photolithography technologies is hinged upon research advances in 
both excimer laser technology and novel materials that possess the required properties 
(high optical transmission at DUV wavelengths, thermal properties, stability after heavy 
DUV radiation exposure) by which the optical elements of the system will be made. 

All illumination systems in optical projection printing tools are designed to provide 

What is known as Knhler illiiminaH'rm r 1-71 Ownl^^tk, • . 

<„ — i L *~j. i-* j K jixwii 6 uiu suwcc 01 an image or me source 

in the front focal plane of the condenser column, the rays originating from each source point 
illuminate the mask as a parallel beam, as seen in Figure 2-2. Each parallel beam is a plane 
wave whose direction of propagation depends on the relative position of the source point 
with respect to the optical axis. Nonuniformity in the brightness of the source points is aver- 
aged out so that every location on the reticle receives the same amount of illumination 
energy. As we will see in subsequent Sections Kohler illumination can be modeled in a con- 
cise mathematical way. 

In addition to dose uniformity, the lithography process should also maintain direc- 
tional uniformity such that the same features are replicated identically regardless of their 
orientations. The shape of the light source is therefore circular (or rotationally symmetric) 
in traditional optical lithography, although this is not true for certain advanced illumination 
schemes such as quadrupole illumination, where directional uniformity is sacrificed in 
order to maximize the resolution of features with certain orientations. 



1 . The term "excimer" originates from the expression excited dimer, meaning a molecule consisting of an 
excited atom and its unexcited (ground) state. Lasing of excimer lasers is highly multimode and this offers a 
crucial advantage for lithography applications, because it relaxes the problem of "speckle" [39] [80] 
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The coherence of the light source is another important attribute. Temporal coher- 
ence 1 is usually not a big concern, since the narrow bandwidth of excimer lasers implies 
high temporal coherence. Spatial coherence (or just plain coherence) on the other hand is 
always carefully engineered and in most cases adjustable. Using special scrambling tech- 
niques, the light emitted from any point of the source is made completely uncorrected 
(incoherent) to the light emitted from every other point. However, light gathers coherence 
as it propagates away from its source [14], [82]. The frequently quoted partial coherence 
factor a is a characteristic of the illumination system and is a measure of the physical extent 
and shape of the light source. The larger the light source, the larger the partial coherence 

factor, and the light source has a lower degree of coherence 2 . In the limit of an infinite 
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Figure 2-2. Illumination with Kohler's method 

Each source point emits a spherical wave that is converted by the illumina- 
tion system into a plane wave incident on the object (photomask). The 
angle of incidence of the plane wave depends on the location of the source 
point (<x,P) with respect to the optical axis (0,0). 



1. Temporal coherence has to do with the ability of light emitted from the same point to interfere with a 
delayed wavefront of itself, whereas spatial coherence has to do with the interference properties of light 
along the same wavefront. 

2. Conceptually, a larger light source contains a greater number of mutually incoherent source points. The 
degree of coherence is therefore less. 
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source, imaging is incoherent and a = ~ . On the other hand, the smaller the light source, 
the smaller the partial coherence factor, and the higher the degree of coherence. Imaging 
with a point source is fully coherent and o=0. Note that a point source in a Kohler illumi- 
nation will result in a single plane wave illuminating the mask and the angle of incidence 
of this wave depends on the relative position of the point source with respect to the optical 
axis. For partial coherence factors between zero and infinity, imaging is partially coherent. 
Typical partial coherence factors in optical lithography range from 0.3 to 0.9. 

2.1.2. Projection system 

The projection system typically consists of a multi-element lens column (up to 30- 
40 lenses) that may also have apertures, filters or other optical elements and it is a marvel 
of engineering precision in order to be able to reliably project images with minimum 
dimensions on the order of lOOnm for state-of-the-art systems. One of the main reasons for 
the required high precision is control of the aberrations, or deviations of the wavefront from 
its ideal shape, but since aberrations are not an immediate topic of this thesis no further 
details on aberrations will be needed. Two relevant parameters of the projection system are 
the numerical aperture, NA, and the reduction factor, R. The numerical aperture is by def- 
inition the sine of the half-angle of the acceptance cone of light-rays as seen from the image 
side of the system. The ratio of image height to object height is by definition the magnifi- 
cation factor M of the system. The inverse of the magnification factor is the reduction 
factor R. Since a typical system in photolithography projects at the image plane a scaled 
down version of the object (mask) 1 , M is less than 1 and R is greater than 1. State-of-the- 
art systems currently have reduction factors of R=4 or 5 and it is exactly because of these 
(relatively) large reduction factors that the powerful methods described in the subsequent 
Chapters have a wide range of applicability. Note that two numerical apertures exist in the 
projection system, namely NA; (or simply NA) and NA 0 , which refer to the half-angle of 
the acceptance cone as seen from the image side and from the object (mask) side respec- 
tively. They are related through the reduction factor as follows: 



ihotomir" 5 haVing 10 d ° Primari ' y re ' aXati0n ° fthC F0CeSS re ^ uirements in «■* fabrication of the 
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Equation 2-1. 



R = 



For a circularly shaped light source the partial coherence factor a mentioned above 
is related to the numerical apertures of both the projection system and the illumination sys- 
tem. Specifically a is given by: 



Equation 2-2. 



a = 



NA C 
NA r 



where NA C is the numerical aperture of the condenser lens (illumination system) and NA p 
is the numerical aperture of the projector lens. Some confusion arises from the fact that in 
the above equation the reduction factor of the imaging system is implicitly taken into 
account. Figure 2-3 clarifies the situation by showing simplified diagrams of two optical 
systems with parameters NA=0.5, a=0.5 and R=5 or R=l. 



5X reduction 



IX reduction 



Illumination optics 
(condenser) 



Imaging optics 
(projector) 




NA=0.5 




Illumination optics 
(condenser) 



NA=0.25 
-29° 



/ NA=0.5 \ 



Imaging optics 
(projector) 



\ NA=0.5 ~/ 

- - - XlUff^J 



Figure 2-3. Numerical apertures and corresponding light 
acceptance cones of illumination and projection lenses for 
NA=0.5, a=0.5 at 5X and IX reduction 
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2.1.3. Photomask 

The photomask, also called reticle, carries the pattern to be printed at a given lithog- 
raphy processing step. The masks of integrated circuits having large die-sizes or footprints, 
(that is, occupying large areas on the semiconductor wafer), typically carry just one copy 
of the chip pattern. A matrix of several chip patterns is contained in one mask whenever the 
chip size permits. Note that the mask is drawn R times the actual size on the semiconductor 
wafer, since the dimensions of the circuit will be scaled down by the reduction factor R. For 
this reason it is not sufficient to just provide feature sizes, since it may not be immediately 
obvious from the context whether these are photomask (object) or resist (image) sizes. A 
typical convention for distinguishing photomask feature sizes from resist feature sizes is to 
include in parenthesis the reduction factor R. For example, a 600nm (4X) line has a size of 
600nm on the mask, and would produce a 600nm/4=150nm line if used in a 4X imaging 
system. Similarly, a 130nm (IX) line refers to the size of a line at the image (wafer) plane 
and would result from the printing of a 130nm line on the mask for a system with R=l, or 
a 520nm (4xl30nm) line on the mask for a system with R=4, or a 1.3um (10xl30nm) line 
on the mask for a system with R=10. 

Depending on their operation principle photomasks can be divided into two broad 
categories: conventional binary or chrome-on-glass (COG) masks and advanced phase- 
shifting masks (PSM). 

A binary or COG mask consists of a transparent substrate (mask blank), covered 
with a thin opaque Him that bears the desired pattern. Light can either pass unobstructed 
through an area not covered by the opaque film or be completely blocked if it is incident on 
an area that is protected by the film. This binary behavior of the transmission characteristic 
of the mask is responsible for its name. The mask blank for DUV lithography typically con- 
sists of fused silica glass that has excellent transmission at *=248nm and somewhat poorer 
but acceptable transmission at X=193nm and A=157nm. The opaque film is typically on the 
order of lOOnm thick and has a chromium (Cr) composition. 

Adding phase modulation to the photomask can profoundly increase the attainable 
resolution. This is the principle followed by phase-shifting masks, which employ discrete 
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transmission and discrete phase modulation . There are many different flavors of PSMs 
depending on the way that the phase modulation is achieved. One of the most promising 
PSM technologies is what is known as alternating phase-shifting mask (alt. PSM, or 
APSM) and since it will be a center part in this thesis it is introduced in some more detail. 
The principle of an alt. PSM is compared with that of a binary mask in Figure 2-4. The 

center line is bordered by transmitting regions with 180° phase difference on an alt. PSM 
and by clear areas of the same phase on a binary mask. The phase difference on the alt. 
PSM leads to destructive interference, resulting in a sharp dark image. The binary mask 

image in not as sharp because of the lack of phase interaction. The 180° phase difference 
is created by etching trenches, also called phase-wells, into the fused silica substrate during 
the alt. PSM fabrication process, which is now more complex than the COG fabrication 
process. The difference in the amount of material removed d etch is such that the path length 
difference between light passing through the different phase regions is half of the wave- 
length in air 2 . The theoretical etch depth difference is determined by the following expres- 
sion: 

n substrate&) d etch ' n air&) d etch = ^ 

A, 

Equation 2-3. d etch = 

2 ^ n substrate^)- "air^)] 

where n substrate (k) is the refractive index of the fused silica substrate at the exposure 
wavelength and n air (k) is the refractive index of air. For commonly used DUV wave- 
lengths (248, 193 and 157nm), n substrat€ (X) = 1.5 and n air (X) = 1 , and Equation 2-3 sim- 
plifies to d etch = X . The performance of alt. PSM depends strongly on the geometrical 

details of the etched phase-wells and accurate characterization and understanding is of 
great importance for the successful utilization of alt. PSMs in the lithographic process. 

1. There are usually two levels of transmission and two, three or four phase levels. Other combinations 
have also been reported. 

2. In principle, any path length difference that is an odd multiple of A/2 would be appropriate, but because 
of fabrication considerations it is almost always A/2. 
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(b) 




Figure 2-4. Comparison of operation principles of an alternating phase shift 
mask (alt. PSM) and a binary (COG) mask 

Cut-planes of geometry of a binary (COG) mask (a) and an alternating phase- 
shift mask (b). The ideal electric field distribution for the binary mask (c) leads 
to a poor image intensity distribution (e) at the image plane, whereas the ideal 
electric field distribution for the alt. PSM (d), because of destructive interfer- 
ence, leads to a robust image. 
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2.2. Resolution in optical lithography 

The smallest resolvable feature or critical dimension (CD) of a lithographic imag- 
ing system is given by: 

Equation 2-4. CD = k x ■ ~ 

where X is the wavelength of the exposure light, NA is the numerical aperture of the pro- 
jector and lq is a process-related factor. In order to increase the resolution (decrease the 
critical dimension) any combination of reducing k 1? decreasing X and increasing NA is 

required. Decreasing X arbitrarily is not possible for many practical reasons 1 , increasing 
the numerical aperture beyond 0.8-0.85 (the physical limit of NA is one) is difficult and 
costly and reducing kj to arbitrarily low values is not possible [103]. Equation 2-4 is some- 
what confusing, in the sense that the ultimate resolution of a system is determined by its 
ability to resolve (small) features that are densely packed and not isolated features. Actu- 
ally, arbitrarily small isolated features can be resolved since in this case there is no theo- 
retical resolution limit [104]! Nevertheless, by Equation 2-4 it is usually implied that the 
distance between the minimum resolvable CD and its neighbors is on the order of one CD, 
that is, the period p of a dense array of (small) features is on the order of two times the crit- 
ical dimension. With that important consideration in mind, theoretically the smallest 
resolvable period p^,, is given by the following expressions for coherent (a=0), partially 

coherent (0 < a < 1 ) and incoherent (a = «>) imaging respectively: 

X 

Equation 2-5. p min = — (a = 0) 



Equation 2-6. p min = 



1 hio<l) 



1+a NA 
2AM' (CT>1) 



1. At the exposure wavelength, a powerful light source needs to exist to satisfy the throughput require- 
ments. Suitable optical materials for lenses and mask substrates and also resist materials are not easy to 
develop at arbitrary wavelengths. 
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Equation 2-7. v . - -— ( a - ^ 

Pm,n 2NA' ( } 

Hence, the critical dimension for 1:1 packed features (feature size of 1CD followed 
by a space of 1CD before the next feature) is half of Pmin and theoretically the minimum kj 
factor is 0.25 (for a>l). However, photolithography imaging and processing below k 1= 0.5 
becomes extremely difficult and it is generally accepted that k 1= 0.3 is a more attainable 
limit. 

2.3. Photolithography simulation and modeling 

The key role of simulation and modeling of the photolithography process is widely 
recognized. In manv asr>ects. simulation in not iust another w»i o.„,~,„_»_ .u 

- J — uiui uugiiiv^nta uic pxoceSS 

development, testing and understanding of the complex relationships of the large number 
of lithographic parameters, but it is a true enabling technology that allows innovative solu- 
tions to be applied. Particularly relevant to the work in this thesis is the role of simulation 
in the application of optical proximity correction (Chapters 3-6). 

One of the distinct modeling phases in photolithography simulation is the simula- 
tion of the photomask imaging. Under specified illumination conditions (size, shape, wave- 
length, bandwidth of source and details of the illuminator), a photomask (binary, alt. PSM 
or other) bearing the pattern to be printed, a given illumination system (NA, filters, aper- 
tures, aberrations) and a certain photoresist "thin film stack" with possible underlying topo- 
graphical features, one is interested in calculating through computer simulation the image 
projected onto (or into) the photoresist. Subsequently, important attributes of the quality of 
this image, such as the critical dimension or the edge slope, can be extracted, or the aerial 
image can be used directly as input to the next module of photolithography simulation that 
models the exposure and development of the photoresist. The same considerations apply in 
the photomask imaging simulation whether it is performed merely to extract the aerial 
image under a certain system setting, or in more complex optimization problems that 
involve, for example, finding the optimum illumination setting that will result in the most 
robust lithographic performance 1 [76], [42], or adjusting the mask pattern such that it will 
result in an image that is a more exact replica of the IC designer's intention 2 . 
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The most rigorous way to simulate the imaging of the photomask would be to 
model the whole photolithography system of Figure 2-1 using Maxwell's equations. The 
continuous form of the two independent Maxwell's equations for linear, isotropic, non- 
magnetic, non-dispersive materials are: 

-> 

* dB 

Equation 2-8. Vx£ = -^ 

-» 

Equation 2-9. Vx// = /+ — 

where J is the electric current density, E, if are the electric and magnetic field strengths 
-> -> 

and D , # are the electric and magnetic flux densities respectively. The following consti- 
tutive relations also apply: 

Equation 2-10. D = e(r)£ 

-> > -> 

Equation 2-11. B = \i(r)H 

Equation 2-12. 7 = a(r)£ 

Every distinct material in the system enters the above equations through its dielec- 
tric permittivity e (or index of refraction n = Je) and magnetic permeability ji, while 
appropriate boundary conditions are enforced at all material interfaces present. Under 

given parameters for the illumination source, the electromagnetic field (E y H) being estab- 
lished everywhere within the system can be determined by the numerical solution of the 
boundary value problem, including the aerial image into or onto the photoresist. However, 
the size of such a problem is monstrous! With a system volume on the order of one cubic 
meter and illumination wavelengths less than \\im the size of the problem expressed in 



1. For example, maximizing the process window or maximizing the overlapping areas of process windows 
for a range of feature sizes. 

2. This process is known as optical proximity compensation (or correction) - OPC. 
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cubic wavelengths exceeds 10 18 ! Clearly, some simplified models need to be considered 
for certain parts of the problem. 

Luckily, both the illumination and the imaging optics parts of the system, as their 
name suggests, can be accurately modeled using the science of optics [11], [36]. Light prop- 
agation through the illumination and imaging optics is modeled with either scalar or vector 
diffraction theory. Elements of the scalar diffraction theory of partially coherent light are 
summarized in the next Section. A treatment of vector diffraction theory can be found in 
[89], [114]. The light propagation effects in the vicinity of the photomask or the resist- 
coated semiconductor wafer may require knowledge of the exact solution of Maxwell's 

equations. Rigorous methods for the calculation of elfi C trnm n <m^ . 

— wutw awancniig mai 

occurs during the passage of light through the object and subject to substrate topography 
during the formation of the image are outlined in Section 2.3.3. 

2.3.1. Simulation of image formation 

The material in this Section follows closely reference [40], but the same ideas and 
concepts can be found in most advanced textbooks on optics. The following theory of 
image formation with partially coherent light will allow numerical calculation of the inten- 
sity distribution expected at the image plane under a specified set of system parameters. 
Integration over the source (Ahhe's method) 

When the illumination of the object originates from a quasi-monochromatic, spa- 
tially incoherent source, as is the case in photolithography systems, there exists a method 
for calculating the image intensity that has the special appeal of conceptual simplicity. 
First, each point on the source is considered individually and the image intensity produced 
by the light from that single point is calculated. Then, the image intensity contributions 
from all points that comprise the source are added, with a weighting factor proportional to 
the source intensity distribution. Simple addition (integration) of the image intensity distri- 
butions is justified, since the original source is assumed to be spatially incoherent. 

Referring back to Figure 2- 1 , under the quasi-monochromatic conditions, each opti- 
cal system can be represented by an amplitude spread function (impulse response). Let 
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F(x,y;a,P) and K(u,v;x,y) be the amplitude spread functions of the illuminating and imag- 
ing systems, respectively. A single source point at coordinates (a,P) emits light that can be 
represented by the time-varying phasor amplitude U s (a,P;t). For Kohler illumination each 
source point corresponds to a plane wave impinging on the object with the angle of inci- 
dence depending on (a,|3), as shown in Figure 2-2. The illumination reaches the object and 
passes through it, resulting in a time-varying phasor amplitude U 0 (x,y;a,f};t) to the right of 

the object given by: 

Equation 2-13. U 0 (x,y;a, P;0 = F(x,y,a, $)T 0 (x 9 y)U s (a, |3;f- 5 X ) 
where is a time delay that depends on (x,y) and (a,|3), and T 0 (x,y) is the amplitude trans- 
mittance of the object, which, for now, is assumed to be independent of the particular 
source point providing the illumination. Finally, the time- varying phasor amplitude of the 
light reaching coordinates (u,v) on the image plane from source point (cc,0) is given by: 

Equation 2-14. 

oo 

tf f (ii,v;a,P;0 = J j K(u, Vy x,y)T o (x,y)F(x,y,0L, &)U s (a, 0;(/- 8, - 8 2 ))dxdy 

-oo 

where §2 is a time delay that depends on (u,v) and (x,y). 

The partial aerial image Ii(u,v;a,P), or the intensity of the light reaching image 
coordinates (u,v) from the source point at (a,f3) is the expected value (time average) of the 
squared amplitude of Uj(u,v;a,P;t). Under the quasi-monochromatic assumption, 
Ij(u,v;a,p) can be calculated to be: 

Equation 2-15. 

CO 

/,(«, v;a,P) = P)J J jjKi^v^y^i^v^y^ 

-oo 

x F(x u y,;a, fi)F*(x 2 , y 2 ;a, $)T 0 (x v yi)T*(x 2 , y 2 )dx x dy x dx 2 dy 2 
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where the asterisk represents complex conjugation and l s (cc,0) is the source intensity at 
(a,3). Finally, the partial intensity IjCu.vja.P) can be integrated over the source coordinates 
(cc,3) giving the result: 



Equation2-16. /,.(„, v) = ///,(«, 0)J J jj K(u , v^, y l)K *(u, v* 2 ,y 2 ) 

-oo -oo 

xF(x vyi ;a, fi)F*(x 2 ,y 2 ;a, ^T^y^T^y^dx.dy.dx^dad^ 

With knowledge of I s , F, K and T 0 the numerical calculation of the image intensity 
distribution is possible. This model of image formation is attributed to Abbe [11 and 
although it is conceptually simple, it is not always the best method to use in practice. 

Representation of the source by an inc id e nt m utual intensity fi mc tion (Hn pli^ mo *u^ } 

Another approach for the calculation of image intensity distributions is possible if 
the explicit integration over the source is suppressed and the effects of the source are rep- 
resented by the mutual intensity function, describing the illumination incident on the 
object. Under the quasi-monochromatic assumption, the time-varying phasor amplitude 
Ui(u,v;t) of the light arriving at image coordinates (u,v) can be represented by the time- 
varying phasor amplitude U c (x,y;t) of the light incident on the object at coordinates (x,y) 
by: 

oo 

Equation 2-17. Ufa**) = j f K(u, v*,y)T 0 (x,y)U c (x,y;t - 8)dxdy 

-oo 

where, again, K is the amplitude spread function of the imaging system, T 0 is the amplitude 
transmittance of the object, and 5 is a time delay that depends on (x,y) and (u,v). Note that 
U Q (x,y;t) and U c (x,y;t) are related through: 

Equation 2-18. U 0 (x,y,t) = U c (x, y;t)T 0 (x, y) 

The intensity at (u,v) is given by: 



Equation 2-19. 



',■<«. v) = <|i/ f (K,v;0| 2 >=> 
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/,•(«, v) = JJJJi«:(ii,vyc 1 ,y 1 )X*(ii > v^,> 2 )r < ,(x 1 ,y 1 )7* <> *(x 2 ,y 2 ) 

-oo 

x (U c (x x ,y x \t- & x )U c *(x 2 ,y 2 \t- & 2 ))dx ] dy l dx 2 dy 2 

Under the quasi-monochromatic assumption the difference |5i-8 2 | is much smaller 
compared to the coherence time 1 and hence: 

Equation 2-20. {U c (x x ,y x \t- S x )U*{x 29 y 2 \t- 8 2 )> = J 0 (x v y { ix 2 ,y 2 ) 
where J 0 is the mutual intensity distribution incident on the object. Finally, the image inten- 
sity is derived upon substitution of Equation 2-20 into Equation 2-19: 

Equation 2-21. 

oo 

'/(".v) = j j !j^^ lt y l )K*(u 9 vv^y 2 )T 0 (x v y 1 )T*(x 2 ,y 2 ) 

-oo 

x J 0 ( x i> y\ y2) dx i d y\ dx 2<ty2 

With knowledge of K, T D and J 0 the image Ij can be calculated. The above nonlinear 
integral equation is often referred to as the Hopkins model and there are situations that this 
model is superior for the numerical simulation of imaging, as will be seen shortly. Note 
that although Equation 2-16 requires six integrations whereas Equation 2-21 only four, the 
latter is not really simpler than the former, for four integrations are in general required to 
determine J 0 . However, when an incoherent source is assumed, as was done for deriving 
Equation 2-16, calculation of J 0 requires only two integrations, for a total of six for the 
complete image calculation. 



1. There exist slightly different definitions of the coherence time x c of a disturbance U(t), involving the 
complex degree of coherence of U(t). Here, it suffices to note that the coherence time is always on the order 
of 1/Av, where Av is the finite bandwidth of U(t). Since U(t) is nearly monochromatic in photolithography 
applications, Av is very small and consequently x c is large. 
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Remarks 

• A key assumption made in the derivation of both Equation 2-16 and Equation 2- 
21 was the independence of the transmission characteristics of the object on the source 
points (O.0). This assumption is not always valid in photolithography systems. While it can 
be easily lifted in Abbe's (Equation 2-16) method it is not trivial to do that in Hopkins' 
method (Equation 2-21), since in the former case the integration over the source is per- 
formed last and a dependence of T 0 on «x,0) is straightforward to include by T o (x,y;a,0), 
whereas in the latter this is not possible because the integration over the source takes place 
first, during the calculation of the incident mutual intensity J 0 . Ways to overcome this sig- 
nificant obstacle are discussed in Section 2.3.5. 

• Under the assumption of a space-invariant or isoplanatic system, which is in gen- 
eral valid for photolithography systems, the following simplifications are true: 

K(u, vjc,y) = K(u - x, v - y), and 

J o( x \>y\*2>yi) = J 0 ( x \ - *2>yi - y 2 ) 

Then, the special cases of completely coherent and completely incoherent illumination 
follow from Equation 2-21. For completely coherent illumination it is 
J o( x i - x 2>y\ - y 2 ) = 1 and Equation 2-21 simplifies to: 

Equation 2-22. /f(l<> v ) = |AT ® T D \ 2 

where ® denotes the two-dimensional convolution operator. The imaging system in this 
case is said to be linear in the complex electromagnetic field amplitude. 

For incoherent illumination it is /„(*, - ^ V] . ^ = 8(X] . ^ ^ . ^ and 
Equation 2-21 reduces to: 

Equation 2-23. /. (m> v) = 0 ^ 2 

and in this case the imaging system is said to be linear in intensity. 
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The Transmission Cross-Coefficient Approach 
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With the assumption of space-in variance Equation 2-21 can be restructured in the 
following significant form: 

Equation 2-24. 

oo 

m s) = J J Tear +/, t + g-j, giT 0 (f +/, + 

The overbars on I x and T Q denote Fourier transforms of the image intensity and the 
object transmission respectively. The quantity TCC is often referred to as the transmission 
cross-coefficient. It is totally independent of the object and is a complete description of the 
optical system from source to image plane. It is given by: 

Equation 2-25. 

oo 

TCCif, g/, g') = J J J 0 (f\ g")F(f+f\ g + g")F*(f+f, g' + g")df'dg" 

-oo 

where, again, the overbars denote Fourier transforms of the respective quantities. 

Finally, the image Ij(u,v) as a function of image-space coordinates can be retrieved 
by the inverse Fourier transform of Equation 2-24. The elegance and power of this 
approach lie exactly in the fact that all illumination and imaging parameters are completely 
modeled in the TCCs, so that they can be pre-calculated and stored via Equation 2-25, and 
the image of a number of different objects can be quickly simulated through Equation 2- 
24. Moreover, the existence of efficient algorithms, like the EFT, for the forward and back- 
ward Fourier transforms involved in this method coupled with methods for the decompo- 
sition of the incident mutual intensity J 0 ([25], [97]) render the Hopkins* method superior 
over Abbe's in a large number of situations. 
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2.3.2. The simulation program SPLAT 1 

A notable implementation of Hopkins' method for the simulation of optical images 
of one-dimensional objects was implemented by O'Toole etal. [61], [60] as part of the pho- 
tolithography process simulator SAMPLE 2 Subsequently, Flanner extended the optical 
imaging simulation program to handle two-dimensional objects [32]. At that time the sim- 
ulator's name was "2D". Toh made significant additions in the simulator's capabilities, 
which included the ability to simulate lens aberrations [94] and also upgraded the name of 
the simulator to SPLAT. Yeung's key theoretical work in extending Hopkins' theory of 
partially coherent imaging to include thin-film interference effects [113] and more accurate 
high-NA models [114], [28] resulted in the current version of the simulation program 
SPLAT, which was used extensively in this work. 

The calculation of images with SPLAT is straightforward. The illumination and the 
projection system are specified through a number of parameters that include X, a, NA, aber- 
rations, reduction factor R, defocus, pupil filters, etc. Next, the photomask pattern is spec- 
ified and this provides the object complex transmission characteristics T Q . Since there is no 
inherent restriction in the complex transmission values, any (periodic) object that adheres 
to a rectilinear grid can be specified, having arbitrary values of transmission and phase. The 
Fourier spectrum To (Fourier series actually, since the object is assumed periodic) of the 
object complex transmission is then numerically evaluated. The TCCs are obtained through 
the numerical evaluation of Equation 2-25 and the intensity at all points in the image field 
is found by taking the inverse Fourier transform of Equation 2-24. 

Other commercial simulation programs with imaging capabilities similar to those 
in SPLAT that are heavily utilized by the lithography research community include PRO- 
L1TH [1 17], SOLID-C [116] and I-Photo 3 . 



1. SPLAT stands for simulation of projection lens aberrations via TCCs. 

2. SAMPLE stands for simulation and modeling of profiles for lithography and etching. 

1 I-Photo is the internally developed photolithography simulator of Intel Corp. It was initially developed bv 
M.chael Yeung m the late 1980's, but the simulator has been through many revisions since then ' * 
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The major limitations of SPLAT (and most other imaging simulators based on 
either scalar or vector diffraction theory) arise from the fact that modeling of the transmis- 
sion of light through the object (photomask) and scattering from substrate topography 
(wafer) at the image plane are circumvented sparing the large computational cost involved 
in properly accounting for these phenomena. The problem of scattering and notching from 
substrate topography has been extensively investigated by other researchers in the past 
[101], [102], [90], [91]. Since it is beyond the scope of this work it will not be considered 
further. However, accurate modeling of the electromagnetic effects that occur during pas- 
sage of light through the object (photomask) is central in this work. Modeling of the object 
in SPLAT is done through what is known as the thin mask approximation, in which the 
photomask complex transmission characteristics are assumed to have ideal transitions at 
all edges of the pattern. In the following, instead of the term thin mask approximation, the 
term Kirchhoff mask model, or k-mask will be adopted, since this approximation stems 
from Kirchhoff s assumed boundary conditions for the problem of diffraction behind a 
screen and is the basis of Kirchhoff 9 s (scalar) diffraction theory. Based on the k-mask 
approximation, when a plane wave is incident on the object, the emerging field is approx- 
imated by a nonphysical piecewise constant modulation. It turns out that when the feature 
sizes on the mask are large compared to the wavelength and when the thickness of the mask 
topography is small compared to the wavelength this approximation is sufficient. How- 
ever, for a large set of objects of interest (phase-shifting masks, masks with OPC) this 
approximation fails and rigorous simulation of the edge effects at the object becomes nec- 
essary. In these cases the tapered on-off transitions of the actual fields at the edges, polar- 
ization dependent effects due to different boundary conditions for the electric and magnetic 
field as well as lateral cross-mixing of the field components in passing through the mask 
cannot be neglected. Simulation has been used to investigate the extent to which such edge 
effects are present and produce significant image effects [105], [33], [2]-[4]. 
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2.3.3. Rigorous simulation of the object and the simulation program TEMPEST 1 

The literature for numerical solution of Maxwell's equations for problems of inter- 
est in photolithography is quite extensive and rich. Research efforts is this area were inten- 
sified after the mid 1980's, as the necessary computing power to handle some problems of 
practical interest was becoming available. The origins of virtually all proposed solutions 
can be traced back to methods that had existed for years, even as early as 1900's, and had 
been applied to the investigation of other electromagnetic phenomena. Such methods 
include Rayleigh's method [73], Green's functions [8], modal expansion methods, the 
method of moments [9] and the finite-difference time-domain (FDTD) method proposed by 
Yee[lll],[92]. 

One of the most frequent classification of rigorous simulations methods in photoli- 
thography is between frequency-domain and time-domain methods. For frequency-domain 
methods the electromagnetic field is expanded into some set of predefined basis functions 
and the unknown coefficients are determined via a matrix inversion problem. The matrix 
inversion part of the problem is computationally intensive and plagues frequency-domain 
techniques, albeit their inherent simplicity. For time-domain methods, a time parameter is 
introduced and the electromagnetic field is found by time-marching, that is, electromag- 
netic interaction with matter is solved in time until the field converges to a harmonic 
steady-state. A comparison of techniques for the simulation of topography scattering can 
be found in [112]. 

The rigorous electromagnetic field solver TEMPEST is used throughout this thesis, 
but it will be obvious from the following Chapters that this work is independent on which 
rigorous simulation method is used for the light transmission through the object. TEM- 
PEST is based on Yee's FDTD method. It was formulated by Guerrieri et al. [44] and 
implemented by Gamelin [34], [35] on a massively parallel computer architecture and it 
was initially used to study the electromagnetic scattering from non-planar topography on 
the wafer (image plane). Reflective notching [91], metrology of polysilicon gate structures 
[90] and alignment mark signal integrity [1 10] were some of the first problems that were 

1. TEMPEST stands for time-domain electromagnetic massively parallel evaluation of scattering from 
topography. - 6 
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tackled with TEMPEST. Subsequently, Wong extended the program to handle three- 
dimensional structures and dispersive materials [105]. He also ported the code to single 
processor architectures and used TEMPEST to investigate issues in photomask edge 
effects, such as the image imbalance in phase-shifting masks [109]. Then, Socha used 
TEMPEST to study imaging and inspection problems and included partial coherence 
effects in scattering from topography [82]. Finally, Pistor optimized and re-parallelized the 
code, added versatile boundary conditions (perfectly matched layer - PML - absorbing 
boundary condition and Fourier boundary condition - FBC) and used TEMPEST to exam- 
ine imaging in extreme-ultra-violet (EUV) lithography and defect printing and inspection 
[68]. 

Using TEMPEST for the simulation of light propagation through the object 

TEMPEST can be used to calculate the time-evolution of the electromagnetic field 
throughout a two or three dimensional structure under the excitation of a monochromatic 
harmonic field. An example 2D simulation is shown in Figure 2-5. A plane wave with 
specified amplitude, phase, polarization and angle of incidence originates at the excitation 
plane and propagates downward into the simulation domain interacting with the mask 
structures that comprise of materials with different (complex, in general) refractive indices. 
The whole transient behavior of the electromagnetic field is readily available as seen for 

example in Figure 2-5(b) that depicts a plot of the real E-field 5 cycles 1 after the excitation 
started. Typically one is more interested in the steady-state electromagnetic behavior, 
hence the simulator is allowed to run for enough cycles until the fields remain constant (to 
within some defined error bound) at the same fraction of each cycle of duration T, in which 
case convergence has occurred. This is seen in Figure 2-5(c) and (d), where the E-field is 
shown after convergence at two time instants that are a quarter-cycle apart. The transmis- 
sion characteristics of the photomask are extracted at the observation plane located below 
the mask. Since the fields are harmonic the phasor notation can be used: 



1. One cycle or period T of the harmonic excitation of frequency f (radial frequency co=27cf) is related with 
the free-space wavelength \$ through: c = X^f = \/T , where c is the speed of light in vacuum. 
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Equation 2-26. E(r, t) = Re^E(r)^ U ] = |E(r)|co>(a» + <j>) 

where E(r) = = |^(?)|cos<I> +>|l(r)| sin0 is the complex (phasor) E-field. If 



the instantaneous field E(r,/)is available at t,=nT (n e N) after convergence and 
t 2 =nT+T/4 (a quarter-cycle apart), then, from Equation 2-26: 

E{r,tJ = |£(r)|cos(ntDr+<|)) = |j?(?)|. C os^ = fle|£(r)J 

t> : \±>.\ f T \ !-> > I 

£tr,r 2 ; = |£(r)|cos^za>r + |+<|>J = -|£(r)|sin<|> = -/m|£(r)J 
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Figure 2-5. Example 2D simulation of the object (photomask) with TEMPEST 

(a) Geometry of a 2D object (mask). A normally incident plane wave originates at the 
excitation plane and propagates downward, interacting with the structure and materi- 
als that are present. The properties of each material are described through the com- 
plex refractive index. The simulation domain is isolated in the vertical direction with 
the PML and is periodic in the lateral direction, (b) Instantaneous (real) electric field 
5 cycles after the birth of the excitation, (c)-(d) Instantaneous (real) electric fields 
after convergence has occurred (15 cycles in this case) that are a quarter-cycle apart, 
(e) Amplitude of the complex electric field. The complex field is constructed using 
the two instant real fields of (c) and (d). (0 Amplitude and (g) phase of the complex 
electric field at the observation plane overlaid with the k-mask model. 



Thus, the complex E-field, E(r) , can be assembled using the following expression: 

Equation 2-27. E(r) = E(r, nT) - jE^r, *T+^j 

The amplitude of the complex E-field throughout the mask of Figure 2-5(a) is 
shown in (e) and the amplitude and phase across the observation plane are shown in (f) and 
(g) respectively. Compare the true transmission characteristics with the ideal k-mask 
model (Kirchhoff approximation for the mask scattering) for this mask, which are overlaid 
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on the plots of (0 and (g). The large error incurred in the subsequent image calculation with 
the k-mask model renders the approximation of the k-mask model insufficient. 

Next, the problem of linking the results of the object scattering acquired through 
rigorous simulation (with TEMPEST) in the equations of the imaging simulator (SPLAT) 
is treated. 

2.3.4. Diffraction orders and integration of rigorous mask simulations in the 
simulation of image formation 

By convention, the object (photomask) is positioned normally to the z-axis and it 
can either be a 2D object (requiring 2D simulation), as in the example of Figure 2-5, or a 
3D object (requiring 3D simulation) that bears a 2D nattem lavont 

a j ~« """"5 in uiv ^ y j-Jiauc. 

Termination of the mask simulation domain in thez-direction is done with the PML absorb- 
ing boundary condition and in the x- and y-directions with periodic boundary conditions 
[71]. The application of periodic boundary conditions to the x- and y- directions of the sim- 
ulation domain has certain implications on the allowed plane wave excitations and the 
resulting angular spectrum representation of the fields across the observation plane The 
electric field (similarly for the magnetic field) of a propagating plane wave is given by (a 

C time-dependence of the steady-state harmonic fields is suppressed): 
Equation 2-28. E (x, y, z) = Eoe^ + ^ + 

where k = (k x , k y , k z ) is the wave vector (units of radians/meter) that indicates the wave- 
length, * = |*| = 2n/X, and Eo is a complex vector indicating the polarization of the 
electric field and also its magnitude and phase. Assuming that the x- and y-dimensions of 
the object are P x and P y respectively, because of the periodic boundary conditions in x and 
y, it is: 



Equation 2-29. E(x + P x , y, z ) = E(x, y, z) => k x P x = m2n => k= m 2* 
and 
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Equation 2-30. £(jc, y + P y9 z) = £(*, y, z) * y P y = n27i fc y = n— 

where m,n are integers. Therefore, periodic boundary conditions restrict the possible 
values of the k- vector of a propagating plane wave inside the simulation domain to a dis- 
crete set, determined by the condition of periodicity of the fields, with periods P x and P y , 
equal to the dimensions of the object in the x- and y-directions. As a result, only a discrete 
set of angles of incidence can exist inside the TEMPEST simulation domain. 

Since only two components of the triad (k x , k y9 k z ) can be independent, one conve- 
nient mapping of plane waves is a k-space diagram (Figure 2-6). With the object periodic- 
ity present in the x- and y- directions, it is convenient to express 1^ with respect to k x and 
k y , which can be arbitrary but adhering to Equation 2-29 and Equation 2-30: 



Equation 2-31. k z = Jk 2 - k 2 x - k 2 y 

Then, the angle of plane wave propagation is given by: 



Equation 2-32. 9 = asin (^j 

and 



Equation 2-33. <|> = atan (^) 

where 9 = Z(k, z) and <|> = Z(k - k • z, x) . The k-space (k x ,k y ) diagram of Figure 2-6 
shows all (discrete number of) plane waves that can exist inside the TEMPEST simulation 
domain. The location of each "x" indicates the direction of propagation through 
Equation 2-31-Equation 2-33. It is implied that the amplitude of the k-vector the 
wavelength) is constant for all points on a k-space diagram. Therefore, points that are fur- 
ther away from the origin correspond to larger angles of propagation with respect to the z- 
axis. Note that the separations of allowable plane waves in k x and k y are in general differ- 
ent, corresponding to different P x and P y . Three important circular limits are indicated on 
the k-space diagram of Figure 2-6: Points within the outermost circle with radius 2tt/A, 
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represent all propagating plane waves. Points within the center circle with radius 
2nNA/kR represent plane waves which get collected by the numerical aperture of the 
imaging optics. Finally, points within the innermost circle with radius 2%oNA/\R repre- 
sent plane waves that are incident on the object. 

t -A k- 
^ x x x x x xx*kx x xx : 

XXX 

xxx: 



xx 




xxxxxxx 



all collected plane waves 
by the imaging system, 
radius 27INA/XR 
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all illuminated plane waves 
by the illumination system, 
radius 2jktNAAR 



Figure 2-6. k-space diagram for plane waves 

The k-space diagram shows plane waves in their spatial frequency representation. The 
origin of the axis represents normal incidence, while points further away from the ori- 
gin represent plane waves with higher angles of incidence (as measured from the opti- 
cal axis). Note that the separation of the discrete number of plane waves allowed by 
the size of the TEMPEST simulation domain can, in general, be different in k x and k 
Also, note that in a typical lithographic system NA/R is < 0.25 so that the radius of 
the circle containing all collected plane waves will be less than 1/4 of the radius of the 
circle of all propagating waves. In this diagram, the relative sizes of these circles are 
drawn different than the typical lithographic imaging situation for illustration pur- 
poses. 



Based on the analysis so far, the following steps can be followed in order to incor- 
porate results obtained from rigorous (TEMPEST) simulations of the object scattering into 
the simulation of image formation: Since the object is illuminated using Kohler's method, 
the radiation from each source point reaches the object as a plane wave, whose angle of 
incidence depends on the location of the source point relative to the optical axis. Although 
the source is continuous, and this implies that an infinite number of plane waves are inci- 
dent on the mask, it can be represented with sufficient accuracy using a discrete rectilinear 
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grid of points. The problem of determining the necessary number of source integration 
points has been worked out by Socha [83]. For each discrete source point a separate TEM- 
PEST simulation can be executed that will determine the true scattering of the object. The 
scattered electromagnetic field across the observation plane is then decomposed into its 
angular spectrum, that is, into a discrete sum of plane waves propagating below the mask 
at angles that accommodate the requirements imposed by the periodicity of the object in 
the x- and y-directions (Equation 2-29-Equation 2-33). This discrete set of complex num- 
bers that fully describe the amplitude and phase of each plane wave component of the elec- 
tromagnetic field scattered by the mask will be hereafter referred to as the diffraction 

orders 1 of the scattered field. A more detailed discussion about plane waves and the angu- 
lar spectrum representation of waves propagating in source-free media can be found in 
Appendix A. 

The next key observation is that the illumination from each individual source point 
is fully coherent. The mutual intensity function is then given by: 

Equation 2-34. 

J 0 ( x \ ' x 2>y\ ' 3> 2 ) = 6 =*J 0 (f>8) = o(f- k x sp ,g- k y sp ) 

where k x sp and k y sp are the x- and y- components of the k- vector of the plane wave inci- 
dent on the mask and they represent the exact location of the source point with respect to 
the optical axis. Therefore the partial image 5lj(u,v) of the respective source point can be 
calculated through Equation 2-24 and Equation 2-25, by taking the inverse Fourier trans- 
form of Equation 2-24. The rigorously calculated diffraction orders are substituted in place 

of the diffraction orders (Fourier transform T 0 ) of the k-mask model in Equation 2-24. 
Finally, the total image Ij(u,v) can be found by repeating the above process for each source 
point and summing up all partial images: 

Equation 2-35. v) = £ 5/.(w, v) 

source - points 



1. Other common names include diffraction harmonics or scattering coefficients. 
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Note that all source points are mutually incoherent, hence all partial images can be 
summed up. This method follows Abbe's formulation of image formation (Equation 2-16). 

Although the aforementioned procedure is simple to implement, it is a total waste 
of computer resources, since it requires a time-consuming (especially for 3D objects) TEM- 
PEST simulation of the mask for each source point of the discretized source and it does not 
take advantage of the fact that the scattering characteristics of the object do not change 
abruptly with the angle of incidence of the incoming plane wave. This problem has been 
successfully addressed in various ways. Wong used a single rigorous (TEMPEST) simula- 
tion to solve the imaging problem of 2D phase shift masks and 3D contact masks [107]. His 
assumption was that the scattering characteristics of the object remain constant for all 
source points for typical imaging situations and can therefore be adequately captured by a 
single rigorous simulation with the plane wave illumination being normally incident at the 
object. This approach was independently verified by Wojcik et al. [100], who concluded 
that for typical lithographic imaging with reduction factors of 4X or 5X there is no need for 
more than one rigorous simulation, even for objects such as phase shift masks that exhibit 
large vertical steps, where the approximation is more suspect. Later, Pi stor re- verified these 
results for binary masks with optical proximity correction and for phase shift masks, and 
showed that this approximation breaks down at a reduction factor of IX that is encountered 
in inspection imaging systems [70]. Alternatively, Socha implemented the Karhunen- 
Loeve [41] expansion as a means of decomposing the illumination of the source into a more 
compact representation [85], rather than the aforementioned plane wave decomposition, 
and Pistor formulated an Abbe-type method [69] that can incorporate as many rigorous 
simulations of different incident plane waves as are necessary to achieve the required accu- 
racy, while avoiding the exhaustive rigorous simulation for all discrete source points. 

For the work in this thesis all imaging simulations that incorporate rigorous mask 
simulations were performed using Wong's approach [108], which is summarized here for 
completeness: 

Using a normally incident plane wave with linearly polarized electric field in either 
x- or y-direction that originates at the excitation plane, the electromagnetic scattering of the 
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photomask geometry is calculated with TEMPEST. The steady-state (after convergence) 
complex electromagnetic field across the observation plane contains all the mask scattering 
information that is necessary. As mentioned before and further explained in Appendix A, 
the Fourier spectrum of each scattered field component (E x , E y , E z , H x , H y , H z ) across the 
observation plane represents a decomposition of that component into a spectrum of plane 
waves propagating at different (discrete) angles below the object. Because of the periodic 
boundary conditions in the x- and y- directions in TEMPEST, the Fourier spectrum is dis- 
crete. The quantity T Q f^ of Equation 2-24 for the k-mask model can be thought of as the 
energy transmitted through the ideal mask. Thus, it is analogous to the Poynting vector 

S = ExH for the rigorously calculated scattered fields. This rigorous mask model will 
be hereafter referred to as the r-mask model. Since by convention the mask lies in the xy- 
plane, the quantity of interest is then the energy travelling in the z-direction. Therefore, 
only the z-component of the Poynting vector S z = (E x H y * - E y H x *) is of interest. The 
Fourier transform of the image intensity can thus be expressed by modifying Equation 2- 
24 as: 

Equation 2-36. 

oo 

m S) = J J TCCif +/, g' + g')[E x (f, g)H*(f, g 1 ) - E y (f, g)H*{f, g')]dfdg' 

-oo 

where the TCC's are given again by Equation 2-25. 

2.3.5. Polarization effects in imaging 

Yet another complication that was suppressed so far arises from the light polariza- 
tion properties. The light being emitted from a laser source in a photolithography imaging 
system can potentially be linearly polarized in a certain orientation, circularly or ellipti- 
cally polarized. As mentioned before, the light from the source goes through a scrambling 
process and as a result each source point emits light that has random polarization with a 
uniform distribution and is also, as seen so far, completely incoherent with all other points. 
Light with random polarization is also called unpolarized or natural light 1 . Viewing the 
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source as an emitter of plane waves that travel in the illumination system, impinge on the 
mask at various angles of incidence, scatter into a spectrum of diffracted orders and subse- 
quently propagate in the projection system, until they reach the image location to form the 
image, is an effective way of modeling and conceptually understanding the whole phenom- 
enon. Since the polarization of each plane wave in the system is random, one needs two 
mutually incoherent, perpendicular directions on the plane normal to k (direction of prop- 
agation) to capture fully the polarization properties of the electric and magnetic fields [13]. 

Traditionally, these directions are the TE ( E 1 (* x z) ) and the TM(Hl(kxz)) defined 
by the following direction (unit) vectors: 

Equation 2-37. = -«in*i+coB#J 

Equation 2-38. a TM = sin9cos<|>i + sin 9 sin <|>y - cosdl 

for a propagating plane wave with k^O (which is the convention used in this work) and 

9 * 0 . When 6 = 0 (normal incidence on the object) a TE and a m are by convention: 

Equation 2-39. ~ ajE = " 

Equation 2-40. a m = x 

The propagation of the waves into the illumination and projection systems typically 
preserves the polarization properties of the waves, although dependencies on wave polar- 
ization could be modeled if needed. The object (photomask) however, more often than not, 
scatters light of different polarizations differently. Moreover, a 3D object illuminated for 
example with a TE wave will scatter energy into diffraction orders of not only the TE polar- 
ization type but also the TM (cross-mixing)! Finally, TE and TM waves recombine in dif- 
ferent ways at the image space to form the image. This phenomenon is more pronounced 
for waves that are travelling at larger angles 9. This is exactly the reason why scalar imag- 
ing theory becomes invalid when modeling the imaging of systems with higher numerical 
apertures, where highly oblique waves are present. 



1. Oftentimes residual polarization exists in the source and it may be necessary to be properly modeled. 
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In summary, the polarization properties are initially set by the light source itself. 
Although the source generally emits unpolarized light, the preferential treatment of the 
illumination/imaging system and/or the object on one type of polarization over the other 
gives rise to polarization effects. The usually unpolarized light emitted by each source 
point can be modeled as two mutually incoherent, linearly polarized plane waves [13] with 
-> . -» - 

E II ajE (TE-polarization) and E II a TM (TM-polarization). Each of these polarizations can 
be treated separately and propagated through the illumination optics that generally pre- 
serve the polarization properties and be incident on the object. A rigorous simulation is 
then executed to determine the scattering properties (diffraction orders) of the object for 
that particular source point and field polarization. Mask polarization effects arise from the 
different scattering response of the object under different incident polarizations. Then, all 
diffraction orders (that generally contain mixed polarization although only TE or TM was 
incident on the object) are propagated through the projection optics, that typically do not 
alter the polarization properties, to the image space where they recombine to form the 
image. Simulation of this recombination process can follow either the simpler scalar dif- 
fraction theory or the more involved vector diffraction theory that properly accounts for 
different field polarizations in the image formation. The two mutually incoherent images 
arising from the two mutually incoherent source polarizations are then averaged to result 
in the final image for an unpolarized light source. This technique is basically followed by 
Pistor [68] in his Abbe's formulation of imaging. Wong [105], [109], Wojcik [100], Socha 
[82] and Adam [2]-[4], [6] used the Hopkins' formulation, calculating the diffraction 
orders for two mutually incoherent, with perpendicular polarizations and normally incident 
plane waves illuminating the object and averaged the two resulting images to calculate the 
total images for unpolarized sources. 

One gray point still remains. Upon close investigation it becomes clear that only 
Abbe's method can properly take into account the source polarization properties. The 

directions of cite anc * cltm change quite abruptly at different source points. This is shown 
in Figure 2-7(a), where the unit direction vectors ajE an d cltm are shown on a k-space dia- 
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gram of the (discrete) source points. Note that a TE always lies on the k x -k y plane, while 
a TM has a z-component out of the k x -k y plane and only its projection on the k x -k y plane is 
shown. Abbe's method can properly account for the rapidly changing directions of the illu- 
minating fields at the various source integration points. On the other hand, in Hopkin's 
method the assumption of constant diffraction orders with respect to the angle of incidence 
of the illuminating wave contradicts the fact that the TE and TM polarizations scatter dif- 
ferently at the object! The TE polarization of a source point lying on the k v =k line for 
example, is seen to have almost equal power in the x- and y-directions. Therefore, even 
within the range of accuracy of Hopkins' approximation, the scattering response of the 
object will be different from source points with k x =k y with TE polarization than the normal 

incidence response for a TE plane wave, where by convention a TE = y ! Here is how this 




all illuminated plane waves 
by (he illumination system, 
radius 2ncrNA//.R 

Figure 2-7. Simulation of the source polarization 

(a) The unpolarized light emitted by each source point can be modeled by two mutu- 
ally incoherent polarizations, a TE and a TMy or any other combination of mutually 
perpendicular directions that are normal to k . The direction of a TE and a TM changes 
abruptly for the various source points and only a TE lies completely in the x-y plane. 
Modeling of the source polarization in Abbe's method is straightforward, but in Hop- 
kin's method the diagram in (a) cannot be applied. Instead, x and y are taken as the 
two mutually incoherent polarizations of each point, as shown in (b). Although non- 
physical, it is a good approximation for typical imaging situations (R=4X or 5X) 
because of the negligible projection of a TM on z . 
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contradiction can be reconciled: In Hopkins' method, TE illumination is modeled by forc- 
ing all source points to emit a polarization that results in a non-physical 1 plane wave inci- 
dent on the object with E y polarization only, same as the polarization for the normally 
incident TE plane wave for which the diffraction orders of the mask have been rigorously 
calculated, as seen in Figure 2-7(b). Similarly, TM illumination is modeled by forcing all 
source points to emit a polarization that results again in a non-physical plane wave incident 
on the object with Ex polarization only, which is the same as the polarization for the nor- 
mally incident TM plane wave for which the diffraction orders of the mask have been rig- 
orously calculated through a second simulation. This is also shown in Figure 2-7(b). The 
validity of this scheme depends on the same conditions that Hopkins' method relies on, Le. 
constant diffraction orders with respect to angle on incidence. In the same way that Hop- 
kins' approximation holds for a relatively small illumination cone, the above decomposi- 
tion between E y and E x polarizations instead of between TE and TM, for all source points, 
is accurate only for source points that are not "too far away" from the origin. The TM illu- 
mination from such source points has negligible E z component and to a good accuracy 
results in a plane wave with E-field completely on the x-y plane of the object space. Then, 
instead of modeling the source polarization by treating each source point as an emitter of 
mutually incoherent TE and TM plane waves, equivalently the mutually incoherent Ey and 
E x polarizations can be used everywhere. Hence, the calculated diffraction orders from two 
separate simulations with the excitation field being a normally incident plane wave with E y 
(TE) or E x (TM) polarization suffice for the application of Hopkins' technique. Again, the 
two aerial images are added for the final image under unpolarized illumination. 

In this work, Hopkins' method with two separate normally incident TE and TM rig- 
orous object simulations is used when calculating aerial images for completely unpolarized 
sources. All approximations are valid for typical imaging situations with R=4X or 5X even 
for NA and a close to unity. The scalar approximations however are not accurate for NA 
exceeding 0.6-0.7, but this is not a limiting factor for this work, since the core of the devel- 



1 . Non-physical because a plane wave propagating in source-free media can never have a field component parallel 
to the direction of propagation k . 
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oped methodology is concerned with speeding up the rigorous simulations of the object and 
can fit equally well in a vector imaging simulator. Extensions of Hopkins' method to cases 
with R=1X or 2X are not obvious. An outline for a possible such extension can be found in 
[113]. The most promising direction for achieving such an extension is believed by the 
author to be a break-up of the source into a number of judiciously chosen, mutually inco- 
herent sectors and the application of Hopkins' method within each sector, where the dif- 
fraction orders of the object should remain relatively constant. Such an extension is beyond 
the scope of this thesis. 

2.4. Approximate vs. rigorous, scalar vs. vector, k-mask vs. r-mask and 
Abbe's vs. Hopkins' methods 

In summary, it is worth stressing the following distinctions about the classifications 
of the various methods that exist in the literature for the modeling of an imaging system 
with applications in projection printing: As mentioned previously, a truly "rigorous" 
method would involve solving the system of Maxwell's equations throughout the imaging 
system. Because of the immensity of such a problem this is not possible. Therefore, strictly 
speaking, all methods that are used in practice for the calculation of optical images are 
"approximate". However, the degrees of approximations used in different methods vary 
and as a result so does the achieved accuracy. Methods that employ a relatively small 
number of approximations are usually advertised as "rigorous" and are typically more accu- 
rate, whereas the ones that involve more approximations are quoted as "approximate" and 
are less accurate. 

The classification of methods between "scalar" and "vector" is a different one. 
Although a "vector" method is in general more accurate than a "scalar" method they are 

both "approximate". In a "vector" method the electromagnetic field (E, H) maintains its 
vector nature throughout the calculations whereas in a "scalar" method at least in one 
module of the calculations the six components of the electromagnetic field (E x , E y , E 2 , H x , 
H y , H z , for a cartesian coordinate system) are treated as scalar quantities and therefore any 
cross-mixing that occurs between those components or any preferential treatment by the 
system on one versus the others is lost. The terms "scalar" and "vector" should be reserved 
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to distinguish whether the illumination and/or the projection system are modeled using 
scalar or vector diffraction theory. 

The terminology of k-mask (short for Kirchhoff mask) and r-mask (short for rigor- 
ously simulated mask) introduced earlier, refers to the way that one models the object 
(mask) scattering. In a k-mask model (also called the thin mask model) no simulation of 
any kind is executed. The mask is assumed to be infinitely thin and the electromagnetic 
field leaving the mask and entering the projection system is approximated by the product 
of the incident field and the ideal complex mask transmission function. In an r-mask model 
(sometimes called the thick mask model) at least one rigorous simulation is needed in order 
to determine the mask scattering. For most imaging situations in practice two rigorous sim- 
ulations are needed to account for an unpolarized source. 

Finally, the classification of methods in "Abbe-type" or "Hopkins-type" has to do 
with the order by which the various integrals in the imaging formulas are evaluated. In an 
"Abbe-type" method, the partial image from each source point is calculated first and then 
all partial images are integrated owing to the assumed source incoherence. This method is 
conceptually attractive and is much easier to incorporate r-mask models, even on imaging 
systems with small reduction factors. In a "Hopkins-type" method a quantity called the 
transmission cross coefficient is calculated first. This quantity characterizes fully all prop- 
agation and image formation phenomena from source to image, except for the mask scat- 
tering effect, that is taken into account separately with either k-mask or r-mask models. The 
key assumption for the application of a "Hopkins-type" method is that the mask scattering 
characteristics do not materially change for the various source points. Owing to this 
assumption the application of the r-mask model in a "Hopkins-type" method for imaging 
situations with small reduction factors, where the assumption of constant scattering char- 
acteristics breaks down, is cumbersome and has not been demonstrated yet. 

The simulation programs utilized (SPLAT, TEMPEST) and properly linked for this 
work result in approximate, scalar, Hopkins' image modeling methods. Both k-mask and 
r-mask models are employed depending on the circumstances of the problem. One of the 
goals of this thesis is to always identify the most efficient simulation strategy for each type 

« 
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of problem and devise methods that can maintain excellent accuracy while allowing for 
large reductions in required simulation time. 
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3 Analysis of OPC Features in Binary 
Masks 



Sub-wavelength lithography, where the size of printed features is smaller than 
the exposure wavelength, places a tremendous burden on the lithographic process. Distor- 
tions of the intended images inevitably arise, primarily because of the nonlinearities of the 
imaging process and the nonlinear response of the photoresist. Two of the most prominent 
types of distortions are the wide variation in the linewidths of identically drawn features in 
dense and isolated environments (dense-iso bias) and the line-end pullback or line-end 
shortening (LES) from drawn positions. The former type of distortion can cause variations 
in circuit timing and yield, whereas the latter can lead to poor current tolerances and higher 
probabilities of electrical failure. 

Optical proximity correction or optical proximity compensation (OPC) is the tech- 
nology used to compensate for these types of distortions. OPC is loosely defined as the pro- 
cedure of compensating (pre-distorting) the mask layout of the critical IC layers for the 
lithographic process distortions to follow. This is done with specialized OPC software. In 
the heart of the OPC software lies a mathematical description of the process distortions. 
This description can either be in the form of simple shape manipulation rules, in which case 
the OPC is referred to as "rule-based OPC" [62], [17]-[22], or a more detailed and intricate 
process model, for a "model-based OPC [54]-[57], [24]-[27], [10], [52], [63]-[64], [74]- 
[75], [87]-[88]. The OPC software automatically changes the mask layout by moving seg- 
ments of line edges and adding extra features that (pre-)compensate the layout for the dis- 
tortions to come. Although after OPC has been performed the mask layout may be quite 
different than the original (before OPC) mask, the net result of this procedure is a printed 
pattern on the wafer that is closest to the IC designer's original intent. 
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The idea of pre-compensation in printing can be traced back to techniques used by 
Eastman Kodak Company [29], who were the first to introduce serifs around corners and 
other "decorative", non-printable, assist features to remove the rounding effects. Sharp 
local reductions of the thickness of the absorbing layer were used in intersecting lines to 
prevent the formation of fillers. These corrections were based on experience gained by a 
process of trial and error. The first attempt for a systematic method for determining such 
corrections appeared in 1981 by Saleh [77], who discussed the general aspect and also 
showed some specific implementations including serifs (corner shape alterations) and sub- 
resolution assist features for better pattern fidelity. Global linewidth biasing for process 
compensation and the use of serifs crafted onto line-ends in memory cells have been 
employed in photolithography for over a decade to enhance the ability to make these struc- 
tures [86]. These techniques were not widely adopted for ICs in general, however, because 
their reliance on (semi-empirical) hand crafted design manipulations made them possible 
only for addressing specific problems in repeating layouts, such as memory cells. 

In the early 1990's the problem of OPC for large mask layouts was addressed for- 
mally as an optimization problem. The size of such a problem is formidable, but through 
introduction of appropriate constraints - induced primarily from mask fabrication con- 
straints, local optimal mask "points" were successfully demonstrated. This research was 
performed independently by Liu, Cobb and Zakhor [54]-[57], [24]-[27], Pati, Wang and 
Kailath [63]-[64] and Stirniman [87]-[88] and has lead to commercially available software 
tools that perform OPC on a full-chip scale. All three approaches rely heavily on speedy 
calculations of the image intensity at selected points of the image field. Although the meth- 
ods by which they achieve this appear seemingly different, they are nevertheless the same, 
in the sense that they all rely on a decomposition of the kernel of the imaging equation for 
partially coherent light. A methodology for improving the insufficient accuracy for OPC on 
alternating phase shifting masks (alt. PSM) is developed in Chapter 4 - Chapter 6. 

The material of this Chapter was originally presented at the SPIE 2000 Microlithog- 
raphy conference. The aim of this Chapter is to develop through simulation a better physi- 
cal understanding of how optical proximity correction works for binary photomasks, and 
how the available geometrical and imaging parameters influence the performance of spe- 
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cific OPC configurations. This is done for "scatter bars" in Section 3.1 and for "corner ser- 
ifs" in Section 3.2. The design approach followed resembles a "rule-based OPC" approach, 
but the developed OPC guidelines could also be helpful for introducing constraints or gen- 
erating good initial guesses in automated "model-based OPC". Then, polarization effects 
of assist features are investigated. OPC shape fidelity and photomask corner rounding are 
also addressed. Finally, the need for rigorous simulations of the photomask scattering (r- 
mask model) as opposed to the simpler k-mask model (in which no mask scattering simu- 
lation is performed) in the presence of wavelength-sized assist features is critically exam- 
ined. 

The analysis methodology used in this Chapter consists of using both k-mask 
models (SPLAT simulation only) and r-mask models (TEMPEST simulation linked to 
SPLAT simulation for the imaging) with two perpendicular polarizations at normal inci- 
dence, as defined in Chapter 2. The k-mask model is quite adequate for binary COG masks. 
However, the level of accuracy is quantified through comparison with the r-mask model. 

3.1. Scatter Bars 

Scatter bars are typically placed near isolated lines as a means of controlling both 
the critical dimension (CD) of the image at best focus, but also improving the through- 
focus robustness of the image [19], [21], [22]. The way they work is shown in Figure 3-1, 
where the resist pattern of the "finger-type" test mask of (a) is depicted in (b) and the aerial 
image calculated with SPLAT in (c). Notice how the isolated part of the middle line and 
also to a lesser extent the outermost lines print fatter in resist. This is not just an artifact of 
the resist, since it is evidenced also in the aerial image in (c). This unacceptable CD varia- 
tion of the middle isolated line can be addressed by using a pair of symmetrically posi- 
tioned lines that are sub-resolution (meaning that the are not wide enough for the imaging 
system to be able to resolve them) and are commonly referred to as scatter bars. The mask 
with OPC scatter bars is shown in Figure 3- 1(d). It is obvious from the resist pattern of 
Figure 3- 1(e) and the aerial image of (f) that the scatter bars significantly improve the 
image of the isolated line, balancing its CD with that of the dense lines. 
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Figure 3-1. Operation principle of scatter bars 

The CD of the center isolated line of mask (a) prints larger than the CD of the dense 
lines as shown in the resist image in (b) and in the aerial image intensity in (c) 
although it is drawn the same on the mask. The addition of a pair of symmetrically 
positioned scatter bars in (d) is shown to significantly improve the resist image in (e) 
and the aerial image in (f). 

The mask designer has essentially three parameters at hand by which she can opti- 
mize the aerial image of an isolated line, namely the size of the scatter bar, the exact loca- 
tion of the scatter bar and the linewidth of the main feature. In this particular study, the 
targeted feature size is 130nm (IX) and the illumination parameters of the optical system 
are: X=193nm, reduction factor R=4, NA=0.7 and partial coherence factor o=0.6, unless 
stated otherwise. All simulations of this Section were performed with SPLAT using the k- 
mask model (no rigorous solution of mask scattering). Figure 3-2 depicts the simulated 
geometry of the isolated line assisted with two symmetrically placed scatter bars. The dis- 
tance from the center of the line to the center of each scatter bar, denoted as d, is expressed 
as the sum of one CD target (=l 30nm) plus pCD^. This implies that if the scatter bars were 
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replaced by lines of width equal to CD target , then p=l corresponds to a 1:1 dense line/space 
pattern (d=2CD target ). 



X=193nm , — H 2 

R=4X Lp {—j Q 

NA=0.7 f ►) |« f 

CD^lSOnm / B CD targe, ^ 



d=(l+p)CD target 



Figure 3-2. Illumination parameters and simulated geometry 



3.1.1. Designing the Placement and Size of the Scatter Bars at Best Focus 

The placement of the scatter bar has a profound effect on the main feature. This is 
shown here by aerial image simulations, in which the only thing that varies is the parameter 
p. The results of these simulations are summarized in Figure 3-3. Observe that when the 
scatter bar approaches the main feature beyond p=l the CD^,^ 1 increases rapidly, until 
the point that the image of the scatter bar completely merges with that of the main feature. 
In the opposite direction, when p»l the line behaves as isolated, as one would expect. For 
intermediate values of p (l<p<3) the placement of the scatter bar controls the CD aerial of 
the line in a fashion captured in Figure 3-3(b), where it is obvious that the CD aerial is more 
sensitive to p for smaller values of a. This behavior appears to be independent of scatter 
bar size, as is also shown in Figure 3-3(b) for SB=55rim and 65nm (IX). 

The second design parameter is the size of the scatter bar. The general rule of thumb 
here is that a larger scatter bar causes a smaller CD aerial of the line for every a (at best 
focus), but for larger values of a the CD aerial becomes progressively insensitive to the scat- 
ter bar sizing. This is depicted in Figure 3-4(a), where the CD aerial vs. SB size is plotted for 
various a. Note that in all these cases the SB is placed at p=l. The slope is steeper for small 
a and flattens out for a > 0.8 . The designer needs to take provision not to oversize the SB, 
otherwise the SB becomes resolvable by the stepper. This fact is depicted in Figure 3-4(b), 

1. CD aeriaI refers to the critical dimension of the aerial image as measured at the 30% level of intensity. This 
corresponds to a simple resist threshold model. 
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Figure 3-3. The effect of the scatter bar placement 

(a) Aerial images at best focus with <r=0.6, SB=55nm(lX) at various placements of the 
SB. (b) Aerial image CD at best focus at the 30% intensity level for various ct and SB 
sizes of 55nm(lX) and 65nm(lX). 

where it is shown that the intensity underneath the scatter bar drops linearly when increas- 
ing the size of the SB. Increasing values of a shift the line upward. Figure 3-4(c) depicts 
the aerial images at best focus for a=0.6, p=l and various sizes of the SB. 

It is also important to understand the role of the partial coherence factor a on the 
effectiveness of the scatter bar. To this extent, a specific design of SB is assumed (place- 
ment at p=l, size of SB=55nm at IX) and a is varied in the simulations. Figure 3-5(a) 
shows aerial images for various o and Figure 3-5(b) shows the changes of CD aerial vs. a. 
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Figure 3-4. The effect of the scatter bar 
size 

(a) Aerial image CD at the 30% 
intensity level vs. scatter bar size 
at SB position p=l, for various a 
(at best focus), (b) Intensity dip 
below the SB vs. SB size at SB 
position p=l, for various a (at 
best focus), (c) Aerial images at 
best focus with a=0.6, p=l for 
various sizes of the SB. 



An important observation here is the fact that the location of the "shadow" of the SB shifts 
with a. For large a it approaches the location of its geometrical shadow. This lateral shift 
implies that in the case that the two scatter bars are misaligned (not symmetrically posi- 
tioned around the feature) the main feature will be displaced! Observe in Figure 3-5(a) that 
the minimum intensity underneath the main feature is increased when a is increased and 
in Figure 3-5(b) that the CD aerial changes significantly with a. For larger a there is a trade- 
off between smaller CD aerial and reduced contrast. 



49 



o 



c 



(a) 



(b) 




: : : : 

— * * i i -: 



— r— -t; 



02 OJ 0.4 01 OS 07 



01 01 03 OA 01 06 07 04 0.9 

Partial Coherence Factor a 



Image position (urn) 

Figure 3-5. The effect of partiai coherence 

(a) Aerial images at best focus for p=l, SB=55nm(lX) for various a. (b) Aerial image 
CD at the 30% intensity level vs. o for SB=55nm(lX) at position p=l (at best focus). 

3.1.2. Out-of-focus Performance of Scatter Bars 

The out-of-focus characterization of scatter bars was also performed with SPLAT, 
k-mask simulations. The results are shown in Figure 3-6. The CD aerial at the 30% intensity 
level with a=0.6 is plotted in Figure 3-6(a) versus defocus for various placements (different 
p values) for SB=55nm(lX). From these plots it is clear that the defocused aerial image is 
poor when no scatter bars are used, but also that the placements at p=1.5 and p=2 do not 
improve the situation. Once the scatter bars are brought nearer to the feature the CD aerial vs. 
defocus curves flatten out and this is indicative of improved through-focus aerial images, 
that resemble the through-focus behavior of dense lines. The best placement in this situa- 
tion occurs at p=0.8. Sizing the scatter bar when p is kept constant does not alter the aerial 
image through focus significantly, as is shown in the CD aerial vs. defocus plots of Figure 3- 
6(b), where p is held constant at 1 and the scatter bar size is increased from 45nm to 75nm 
in steps of lOnm. 

In summary, the design process and verification of performance of scatter bars next 
to isolated features could proceed as follows: First, according to the illumination conditions 
(X, NA, a) the scatter bar size is chosen from a design graph similar to that of Figure 3-4(a) 
that is consistent with the discretization of the mask making tool and results in the required 
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Figure 3-6. Through-focus characterization of scatter bars 

Aerial image CD at the 30% intensity level vs. defocus at a=0.6 for (a) 
SB=55nm(lX) at various positions p and (b) position p=l for various sizes of the 
SB. 

correction of the CD aeria i. Care must be exercised, so that the scatter bar is below the 
camera resolution and will not print (using Figure 3-4(b)). Second, from Figure 3-3(b) the 
acceptable range of p is determined so that the CD aeria i remains close to its minimum for 

p = 1 , say 0.75<p<1.2. Third, through aerial image simulation the position p (0.75<p<1.2) 
is found that shows the best performance through focus. Figure 3-6(b) gives convincing 
evidence that the size of the scatter bar does not affect the out-of-focus performance. Any 
residual difference from the target CD aerial will have to be corrected by biasing the main 
feature. (This makes things a bit more complicated, because biasing the feature changes p 
and iteration will be needed). Finally, the ability to adjust the partial coherence factor a 
provides one extra knob to control the CD aerial , as is apparent from Figure 3-5(b). 

3.1.3. Polarization Effects of Scatter Bars 

The simulations in the above discussion were performed with SPLAT, using the k- 
mask (Kirchhoff mask) model, which is based on "vertical propagation", where the mask 
is assumed to vertically transmit the jncident field. It is of interest to know, how well the 
above assumptions agree with the rigorous, r-mask model, where the mask scattering is rig- 
orously simulated with TEMPEST, and also how scatter bars near isolated features behave 
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at different polarizations of the incident light. Figure 3-7 depicts a sample TEMPEST sim- 
ulation where p=l and SB=50nm(lX). The near field amplitude at a zx cut-plane is shown 
in (a) for TE polarized incident light (electric field parallel to the line) and in (b) for TM 
polarized light (electric field perpendicular to the line). Observe that in the TM case the 
electromagnetic field underneath the mask corresponds to a travelling scattered wave, 
whereas in the TE case the field also exhibits a standing wave pattern along the x-axis. This 
can be explained as follows: In the TE case there exist induced currents onto the chromium 
layers, that are aligned with the incident radiation. These currents at the edges of the feature 
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Figure 3-7. Amplitude of the complex near field 

(a) TE polarization - E y only and (b) TM polarization - E x only. In this example p=l 
and SB=50nm(lX). Note that because of the symmetry the TEMPEST simulation 
is performed only on the right half of the mask geometry 
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and the scatter bars re-radiate cylindrical waves, which interfere constructively or destruc- 
tively, and result in the peaks and valleys along x that are spaced by ~100nm {i.e. ~\I2). In 
the TM case the induced currents are weaker, because the field is now perpendicular to the 
line and they cannot redirect energy in the region underneath the mask. This results in the 
transmitted field being basically just the incident field in the open areas. 

The aerial images for the TE and TM cases are compared for this example simula- 
tion in Figure 3-8, together with the SPLAT simulation. The aerial images when no scatter 
bars are used are also shown for comparison. Observe that the scatter bars and the main 
feature appear wider in TM excitation than in TE and narrower with SPLAT simulation. 
Figure 3-8(b) depicts the intensity dip underneath the scatter bar for various sizes of scatter 
bars for the TE, TM and SPLAT cases. This dip is always larger (smaller intensity) for the 
TM case and smaller (higher intensity) for the SPLAT case and the TE case falls in 
between. Observe also that for larger sizes of scatter bars the three cases give even more 
similar results. For smaller sizes of scatter bars a simple perturbational model can be used, 
in which the intensity dip is approximated by the square root of the intensity dip of an 
opaque pin spot. Although there is a difference of up to 10% between the TEMPEST and 



(a) (b) 




Image position (^m) Size of SB (X/NA) 

Figure 3-8. The effect of polarization on the performance of scatter bars 

(a) Aerial images at best focus without OPC and with a SB=50nm(lX) at p=l for 
TE and TM polarizations. The SPLAT simulation result is also shown. 
(*=193nm, NA=0.7, a=0.6, Mag.=4X, CD target =i30nm(lX)) (b) Intensity dip of 
the scatter bar for various scatter bar sizes as predicted by SPLAT, TEMPEST 
(TE and TM polarizations) and a simple perturbational argument, when SB is 
placed at p=l. 
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SPLAT simulations, the general guidelines for the proper design of the scatter bars dis- 
cussed in Sections 3. 1 . 1 and 3. 1 .2 still hold, but a final fine-tuning using the rigorous TEM- 
PEST simulation may be necessary, especially when there is need to account for 
polarization effects. 

3.2. Serifs 

Square serifs are typically introduced at the corners of lines in order to locally cor- 
rect the aerial image and in an effort to reduce the line-end shortening (LES). An example 
of a mask pattern before and after the application of OPC serifs at interior and exterior cor- 
ners is shown in Figure 3-9(a) and (b). The aerial image is seen to significantly improve in 
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Figure 3-9. Operation principle of corner serifs 

Mask before (a) and after (b) corner serif OPC. The aerial images are shown in (c) 
and (d) respectively and the 30% intensity contours are overlaid with the ideal (geo- 
metrical) image. The shape fidelity near the corners is much better after the OPC and 
the line-end shortening (LES) is significantly reduced. 
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3.2. 1. Modeling OPC serifs through simulation 

In this Section the effect of the size and placement of the OPC serif on the LES aerial 
correction of the aerial image is modeled. Then, the additional accuracy from rigorous 3D 
electromagnetic simulation (TEMPEST) in predicting the correction in the LES aeria] intro- 
duced by a specific OPC design is determined. To do this, a number of cases that are sim- 
ulated with both TEMPEST and SPLAT are compared. In the following, the targeted 
feature size is 130nm (IX) and the illumination parameters of the optical system are: 
X=193nm, reduction factor R=4, NA=0.7 and partial coherence factor a=0.6. Since only 
the aerial image improvements near the end of the isolated line are of interest, no scatter 
bar OPC is used for linewidth correction and out-of-focus improvement of the aerial image. 
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Figure 3-10. Different OPC placements for a clear field (CF) mask 

(a) diagonal non-overlapping, (b) 3/4 diagonal overlapping, (c) side and (d) top. 
Note that because of the 2-fold symmetry with respect to the x, y axes only the top 
right quadrant of the simulation domain is shown. 
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Figure 3-10 depicts four different OPC placements for a clear field (CF) mask: (a) 
diagonal non-overlapping, (b) 3/4 diagonal overlapping, (c) side and (d) top. Intermediate 
cases where the OPC is shifted along the diagonal to render different than 3/4 diagonal 
overlapping or up-down for different side and left-right for different top placements were 
also simulated. In Figure 3-11 an example 3/4 diagonal overlapping design is compared 
with the reference case (no OPC used). The side of the square serif is 65nm(lX) and it is 
positioned so that the center of the serif is located at the comer of the line. This geometry 
is simulated with TEMPEST (r-mask) and SPLAT (k-mask) and the contours of the 30% 
clear field intensity of the aerial images are shown in Figure 3-1 1(b). In the same plot the 
30% contours for the reference case are also drawn, which again is simulated with TEM- 
PEST and SPLAT. The ideal rectangular image is also shown for comparison. Note that 
both simulators predict a ~10nm linewidth bias on each side of the CD aerial , which is the 
approximately the same value with that resulting from Figure 3-3(b) for p=5 (nearly iso- 
lated line) and o=0.6, but they differ in their line-end predictions. The LES aeria] resulting 
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Figure 3-1 1. An example of 3/4 diagonal overlapping OPC design is compared with the 

reference case (no OPC used) 

(a) The side of the square serif is 65nm(lX) and it is positioned so that the center of the 
serif is located at the corner of the line, (b) Contours of the 30% clear field intensity of 
the aerial images from TEMPEST and SPLAT simulations. The ideal rectangular 
image is also shown for comparison. 
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from a TEMPEST (r-mask) simulation is always smaller (i.e. longer line) than that result- 
ing from a SPLAT (k-mask) simulation of the same CF (clear field) geometry. But what is 
a more meaningful comparison, is rather one between the correction in the LES aerial from 
the reference to each OPC case resulting from the TEMPEST and SPLAT simulations of 
the same CF masks. In the example shown in Figure 3-11 this correction comes out to be 
31nm or 34nm from the SPLAT or TEMPEST simulations respectively. 

In Figure 3-12 the square root of the correction of the LES aerial is plotted vs. the 
equivalent size of OPC, where the equivalent size is the square root of the non-overlapping 
area, for all the different OPC serif sizes and placements that were simulated with TEM- 
PEST and SPLAT, and for both mask polarities (dark field, clear field). For a CF mask the 
correction in the LES aer j a i predicted by TEMPEST is larger by -9% than that predicted by 
SPLAT, for all OPC designs, whereas for a DF mask TEMPEST predicts a smaller 
LES aerial correction by -6% compared to SPLAT, for all OPC designs. Why this is hap- 
pening can be explained by invoking the boundary conditions of the electromagnetic field 
in the vicinity of the serif: The electric field that is parallel to the square serif sides (as is 
true for good conductors) has to vanish within the chromium, which requires the tangential 



equiv. size of OPC in A/NA (1 X) 

Figure 3-12. Square root of the correction of the LES aerial vs. equivalent size of OPC 

The equivalent size is the square root of the non-overlapping area for all the different 
OPC serif designs considered. Both mask polarities (dark field, clear field) are exam- 
ined. The points that do not fall on the respective straight line correspond to OPC 
designs where the serif was placed in highly off-diagonal locations. 
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fields in the neighboring open areas to rapidly vanish in close proximity to the chromium. 
This effectively implies that the electromagnetic size of a dark serif (CF mask) is somewhat 
larger than its geometrical size and that the electromagnetic size of an open serif (DF mask) 
is somewhat smaller than its geometrical size. Hence a dark serif will actually be responsi- 
ble for a bigger correction in the LES aerial (ALES aerialTEMpEST >ALES aerial)SPLAT)> 
whereas a clear serif for a smaller correction (ALES aeria] )TEMPESX <ALES aerial SPLAT ). The 
above trend with serif area can be further pursued by using perturbational arguments to 
relate the size of a small pin-hole with the peak intensity of light that goes through, and also 
the size of a small opaque pin-spot with the minimum intensity that appears at its shadow. 
In that respect, from a previously proposed model [59]: 

Equation 3-1. / = s 5^ CD Y 

peak- PH ■ yX/(NA)J 

and using perturbational arguments: 

Equation 3-2. . , . 2 JT^ 

For example the 65nm x 65nm (IX) square opaque pin-spot used in the simulations 
shown in Figure 3-11, when simulated in an isolated configuration with TEMPEST (r- 
mask), results in an intensity dip of 0.64, whereas the intensity dip resulting from the 
SPLAT simulation is only 0.69. This renders an equivalent electromagnetic size of 69nm x 
69nm from the r-mask model or 63.5nm x 63.5nm from the k-mask model, i.e. a difference 
in their side of -8%, which correlates with the -9% difference in the effectiveness of this 
serif when placed near the corners of the isolated line. An open square pin-hole of the same 
dimensions appears smaller by -5.5% in TEMPEST than in SPLAT. Therefore, its effec- 
tiveness is reduced when it is placed near the comers of an isolated open line (DF mask), 
and the reduction in the correction of the LES^^ is -6%, as shown in Figure 3-12. 

Note that in Figure 3-12 the points that do not fall on the respective straight line cor- 
respond to OPC designs where the serif was placed in highly off-diagonal locations, but 
even in those situations the r-mask to k-mask (TEMPEST-SPLAT) correction bias is close 
to- the 9% (CF) and 6% (DF) values that were found above. The important result that the 
graph of Figure 3-12 conveys is exactly the fact that a rigorous electromagnetic simulation 
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is not necessary in order to predict the correction in the LES aeria] that a certain OPC serif 
design will cause. This can be done with the k-mask model (just SPLAT) just as well, but 
one needs to account for the 9% (CF), or the 6% (DF) r-mask to k-mask difference. Viewed 
a little differently, in order to correct for a certain amount of LES aeria j, just SPLAT simu- 
lations can be used for the design process, but in the case of a CF mask a correction for 9% 
smaller LES aerial needs to be done and for the case of a DF mask a correction for 6% larger 
LES aeria j. 

3.2.2. Corner Rounding 

Here, isolated lines are considered on both clear field and dark field masks with 
rounded corners. The reason why corner rounding of the mask feature is examined, is 
because its aerial image behavior with respect to the degree of roundness resembles the 
OPC serif design, where now the roundness of the comer can be thought of as an anti-serif, 
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Figure 3-13. The effect of corner rounding 

Amplitude of the near field of a reference (not rounded) and a rounded 
case (a) for a CF mask and (b) for a DF mask, (c) Increase of the 
LES aeria j vs. radius of curvature for CF and DF masks. 
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as it will be obvious shortly. Figure 3-13(a) depicts the amplitude of the near field of a ref- 
erence (not rounded) and a rounded case for a clear field mask and Figure 3-13(b) for a dark 
field mask. In both (a) and (b) the incident light is parallel to the line (TE). Varying the 
radius of curvature of the rounded mask feature the resulting aerial image is examined and 
the amount by which the contour of 30% intensity retreats from the respective contour of 
the ideal (right angled) mask is deduced, that is the increase of the LES aerial . The increase 
of the LES aeria , with respect to the reference is plotted versus the radius of curvature for 
both clear and dark field masks in Figure 3- 13(c). From these plots it is apparent that the 
LES aerial increases quadratically versus the radius of curvature, which implies that it is pro- 
portional to the area missing from the corner because of the roundness, hence it is valid to 
think of comer rounding as anti-serif OPC in a CF mask or serif OPC in a DF mask, which 
lead to the opposite than the desired direction of LES aeriaJ reduction. 

3.2.3. Shape fidelity requirements for good performance of OPC serifs 

Knowledge of the required shape fidelity for the OPC serif is of paramount impor- 
tance to the mask designer and especially to the mask manufacturer. From the above dis- 
cussion it is apparent that the increases or decreases of the LES aerial correlate linearly with 
the feature area that is missing or is extra near the corners. If that is the case, why should 
the mask manufacturer invest in time and effort to faithfully replicate "perfect" square OPC 
features onto the mask? Is it possible to get by with an easier to print OPC serif shape that 
still introduces the same non-overlapping area? The answer to the above questions, as far 
as aerial image simulation predicts, is that indeed, it is the extra area introduced by the serif 
that is responsible for the correction in the LES aeria] , whereas the details of the shape are 
second order important. Figure 3-14 depicts one example, where a square OPC of side 
0.1A7NA is compared with a "mouse ear" (circular) OPC of radius 0.06A/NA. They are both 
placed with their geometrical centers coinciding with the comers of the isolated line, there- 
fore having equal non-overlapping areas. Figure 3- 14(a) depicts the amplitude of the near 
field for the reference case, the square OPC and the "mouse ear" OPC, (b) depicts the aerial 
image intensity contours for the two OPC shapes and in (c) a cut-line along y at x=0 is 
shown, from which it is obvious that the two aerial images are nearly identical. In Figure 3- 
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15 the LES aerial correction the non-overlapping area of the serif is shown for DF masks and 
various sizes of square and "mouse ear" serifs, when their center is aligned with the corner 
of the main feature. 

(a) 

Square OPC 



|EJ Near Fields 



Reference 
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"Mouse Ear" OPC 
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Figure 3-14. Effectiveness of different OPC shapes 

(a) Amplitude of the near field of the reference (no OPC), a square OPC with side 
0.1A7NA and a "mouse ear" OPC with radius 0.06A/NA, so that they have the same 
non-overlapping area, (b) Aerial image intensity contours of the square and "mouse 
ear" OPC. (c) Cut-line of the aerial images along the y-axis at x=0. 
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Figure 3-15. Requirements for OPC shape fidelity 

LES aeria | correction vs. the non-overlapping area of the serif for DF masks 
and various sizes of square and "mouse ear" serifs, when their center is 
aligned with the corner of the main feature. 
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3.3. Conclusions 

Scalar imaging and rigorous electromagnetic simulation have been used to generate 
quantitative design data on the placement and sizing of OPC features. Studies of scatter 
bars showed that both the best-focus and out-of-focus aerial image CD are more sensitive 
to the scatter bar placement rather than size. For very high c the CD aerial is a weaker func- 
tion of both placement and size of the scatter bars. The r-mask to k-mask (TEMPEST- 
SPLAT) difference for scatter bar simulations is typically only -10%, but the same trends 
for scatter bar placement and size still apply. 

Using a perturbation approach, it is shown that OPC serifs can be adequately char- 
acterized by SPLAT (k-mask) simulations onlv. The r-masV tn v.™^ M .„u.> u: w 

sponds well with the bias of the two models in the case of small isolated holes and posts. 
This behavior is accurately characterized for diagonal placements of the serifs, but for off- 
diagonal placement it is not as accurate. The actual shape of the serif is shown by rigorous 
3D electromagnetic simulation to be less important than previously expected. It is shown 
that the correction of the LES aerial is proportional to the non-overlapping area of the serif. 
Corner rounding can be regarded as an anti-serif case in clear field masks and serif in dark 
field masks. The amount of LES is proportional to the area missing because of the rounding 
of the mask feature. 
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^/fl Domain decomposition methods for 
■ the rapid simulation of photomask 
scattering 



SIMULATION has been very successful in the design and optimization of phase- 
shifting masks. The infamous intensity imbalance problem of an alt. PSM was first discov- 
ered by Wong using simulation [105], [109] and later verified experimentally [66], The 
optimization of the trench design of alt. PSMs through simulation was demonstrated by 
Friedrich et al [33]. Hotta et al included the spatial frequency spectrum view in the anal- 
ysis of an alt. PSM [46]. Hotta also revealed the role of the Chromium-based absorption 
layer that was first noted by Yuan [115]. One of the remaining major challenges, that is 
addressed in this Chapter and further developed in Chapters 5 and 6, is the 3D simulation 
of large non-periodic mask layouts and the considerable computational effort that is 
involved. The challenge is to not only get the near fields, but to obtain an equivalent and 
compact representation of the mask that fully captures its diffraction properties. If such a 
task can be completed rapidly, then it can play a major role as a part of a model-based OPC 
method that can accurately handle alt. PSM (or any other advanced photomask). The use 
of an OPC software tool with this important capability can modify the mask in such a way 
that the resulting spectrum meets the design goals. Issues such as pitch dependent focus 

behavior, forbidden pitches and small MEEF 1 in alt. PSMs [81] could thus be addressed 



1 . MEEF stands for mask error enhancement factor. It is typically defined as: 

MEEF - ACD wafer /ACD mask , where ACD wafer is the variation of the critical dimension of the (pre- or 
post-etch) printed image on the wafer from a nominal CD wafer and ACD^^ is the variation of the critical 
dimension on the mask from a nominal CD,,^. The nonlinearities of the process lead to a widely varying 
MEEF at various nominal CDs and this is a serious problem. An intriguing look at the MEEF as a matrix 
that allows two-dimensional proximity effects to be quantified can be found in [43]. 
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effectively at the mask simulation level, without rigorous techniques that can only be used 
sparingly. 

This Chapter progresses as follows: First, building upon an example that illustrates 
the discrepancy between the two imaging simulation approaches, namely the k-mask model 
and the r-mask model for alt. PSMs, the new quasi-rigorous mask model (qr-mask) that 
relies on a domain decomposition method (DDM) is introduced. The linearity rooted into 
the Kirchhoff-Fresnel diffraction integral is used as a starting point for the domain decom- 
position method. The qr-mask model is further developed into the more compact matched- 
Kirchhojf mask model (henceforth mk-mask), that allows for the rapid and accurate simu- 
lation of the imaging problem. The extension of the qr-mask and mk-mask models in three 
dimensions is demonstrated and benchmarking results that demonstrate the validity of the 
method for 3D simulations are presented. 

Two other important issues are also discussed at some length, namely the applica- 
bility of this method for arbitrary polygonal shapes of the Manhattan-type that are generally 
encountered in mask layouts and the extension of the method in imaging situations where 
the mask diffraction orders are not constant for all angles of incidence within the illumina- 
tion cone. The important topic of cross-talk between adjacent features in alt. PSMs is an 
integral part for the general applicability of the models developed in this Chapter and is 
covered in the next Chapter. 

4.1. Background 

4.1.1. Superposition and domain decomposition 

Inherent in the linearity of Maxwell's equations for the electromagnetic field is the 
principle of superposition. One needs to be careful though, in the fact that the superposition 
principle applies to the sources of the electromagnetic field and not to the geometry or the 

materials present in the problem. For example, if the tangential En and H t \ fields are 
given everywhere along a closed surface that surrounds a volume of materials e u a 
(where the material properties are a function of position), then a unique solution exists 
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everywhere within the volume, say E a \ and Hal. Now, if a different set of tangential 

fields, Eti and ///2 , are enforced (owing to different external sources) along the same sur- 
face enclosing the same geometry of materials e fl , \i a , o a , a new unique solution is being 

-> 

established everywhere within the volume, say E a 2 and H a 2 . By the linearity of Max- 

-> -» 

well's equations the superposition principle follows, according to which, if E t \ + Eti and 
-» -> 

#m +Ha are the tangential fields enforced along the surface, the solution within the 

■> -> -> -» 

volume will now be the sum of the two solutions, namely E a \ + E a 2 and H a \ + #a2 . How- 

(a) 






£rl + E,2, Ht\+Ht2 



(b) 






Figure 4-1. Superposition principle 

The superposition principle follows from the linearity of Maxwell's equations 
with respect to the sources of the EM field, as shown in (a). However, a superpo- 
sition of materials (geometries) such as the one shown in (b) is false. 



65 



t c 

ever, if the tangential En and Ml fields are given everywhere along the same closed sur- 
face that surrounds a volume of new materials e b , \i b , a b (different geometry), then a new 

■* -* 

unique solution exists everywhere within the volume, say Eb] and Hb\ . Note that the 

•* -» 

superposition principle does not apply in this case, that is, if Et\ and Ha are the tangential 

fields enforced along the same surface that surrounds a volume of materials 

E a + E b> V a + °a + a b ( com P osit e geometry), it is easy to prove that the solution within 

-» -» -» -» 

the volume is not £ a i + Ebi and # a i + Hb\ , as one might be tempted to suggest. 

A different kind of superposition will be invoked in the following, namely the 
superposition stemming from the linearity of the Kirchhoff-Fresnel diffraction integral 
[15], that gives the (scalar) disturbance U (any field component) at a point P behind an infi- 
nite screen bearing an opening W, when illuminated from a distant point source at P 0 
located on the other side of the screen: 

Equation 4-1. U(P) = J J^(l + cosx)</S 

0 w 

where X is the wavelength, A is a constant (complex amplitude from the source), k is the 
wavenumber and all other parameters are depicted in Figure 4-2. The diffraction integral of 
Equation 4-1 extends over W, the screen's opening. If the screen consists of multiple open- 
ings Wj, according to Equation 4-1, one could calculate the disturbance Uj(P) due to each 
opening separately and then sum the partial disturbances for the total U(P): 




/ 

/ 



Figure 4-2. Illustrating the Kirchhoff-Fresnel diffraction formula 
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Equation 4-2. 



U(P) = £t/,(P) 



The above equation is the essence of the domain decomposition that will be applied 
in the following in alt. PSMs. Every single opening of an alt. PSM will be considered (sim- 
ulated) separately and the scattered field below the mask will be found via the DDM as the 
sum of the scattered fields below each opening. This is shown in Figure 4-3, where an arbi- 
trary alt. PSM layout (2D or ID) is decomposed via the DDM into its constituent single 
opening masks. For now, each 2D opening is restricted to be rectangular. Arbitrary polyg- 
onal shapes are addressed in Section 4.3.5. Note that according to the simulation strategy 
developed in Chapter 2, only the scattered field from a single, normally incident plane 
wave needs to be calculated. 
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Figure 4-3. The domain decomposition method (DDM) applied in 2D layouts 
(3D masks) and ID layouts (2D masks) 
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4.1.2. Simulation approaches and the imbalance problem in alt. PSM 

The simulation of an image intensity distribution (aerial image) that is projected by 
an optical system through a photomask onto (and into) the resist-coated semiconductor 
wafer can be performed following one of the two distinct paths shown in Figure 4-4, as dis- 
cussed in Chapter 2. In the more computationally laborious path of Figure 4-4(a) the exact 
three-dimensional geometry of the photomask is used for a rigorous calculation of the elec- 
tromagnetic field (EMF) that is being established as a result of a certain illumination. Such 
a rigorous solution demands large computational resources, since even for a minute 10X x 
10A layout more than 100Mb of memory and lOhrs of simulation time may be required. 
Obviously, large, arbitrary mask layouts cannot be handled. The complex (amplitude and 
phase) steady-state near EMF that is scattered below the mask constitutes the r-mask 
model, which can be thought of as a rigorous complex mask transmission function (CMTF) 
that, owing to the field component cross-mixing, consists in general of all 6 field compo- 




EMF (r-mask) 




CMTF (k-mask) 



diffraction 
spectrum 



a*- projection lens 




aerial image 



_A_ 



Figure 4-4. The two distinct simulation approaches for the aerial image 
calculation of an alt. PSM 

In (a) the near electro-magnetic field (EMF) is found through rigorous simu- 
lation of the exact mask structure (with TEMPEST) and used as the r-mask 
model for the aerial image calculation (with SPLAT), while in (b) a complex 
mask transmission function (CMTF) is assumed and used as the k-mask 
model for the aerial image calculation (with SPLAT). 
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ncnts (E x , E y , E z , H x , H y , H z ), that is, one CMTF per field component. As seen in Chapter 
2, each complex field component across the observation plane can be decomposed into its 
complex angular spectrum (Fourier spectrum of the near EMF). Note that this diffraction 
spectrum is continuous for a completely isolated mask structure and discrete for periodic 
mask structures. Typically, most simulators use periodic boundary conditions in the lateral 
mask directions and subsequently this makes the spectrum discrete. The frequency separa- 
tion of the complex discrete diffraction orders depends solely upon the mask period in x- 
and y-directions, as seen in Chapter 2. 

In photolithographic imaging for IC, the objective lens acts as a low-pass filter and 
only the center part of the diffraction spectrum is collected and contributes in the aerial 
image formation. In the case of a discrete spectrum, the diffraction orders that are collected 
by the imaging system's acceptance cone (NA) are referred to as through- the -lens (TTL) 
diffraction orders. Specifically, the collected plane waves have wavenumbers such that: 

Equation 4-3. Jk x + k y < y ^ 1 

where X is the illumination wavelength, NA is the numerical aperture of the projection 
optic, a is the partial coherence factor of the illumination and R is the reduction factor of 
the optical system. 

In the faster simulation approach of Figure 4-4(b) the rigorous but time-consuming 
calculation of the photomask EMF is omitted and, instead, a scalar CMTF is assumed, i.e. 
a k-mask model is used that typically consists of a piecewise constant function coinciding 
with the photomask transmission properties in x, y. The limits of validity of this approach 
have been investigated extensively by many authors and it has been established that it 
breaks down when either the minimum feature on the photomask layout or the vertical 
structure of the photomask is comparable to the illumination wavelength. Unfortunately, 
this is true for advanced photomasks with OPC or alternating PSM and the simulation path 
of Figure 4-4(b) leads to inaccurate results. It suffices to say that the accuracy of this 
approach is contingent upon how closely the spectrum resulting from the k-mask model 
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matches the spectrum from the rigorously calculated r-mask model within the bandwidth 
of the projection system, as depicted pictorially in Figure 4-4(b). 

An example of a CD=400nm (4X), 1:1 line/space, 0°/180° alt. PSM is shown in 
Figure 4-5(a) and the resulting EMF solution calculated with TEMPEST for X=193nm TE- 
excitation 1 is shown in Figure 4-5(b). The amplitude and phase of the near EMF at the 
observation plane located 50nm below the end of the 80nm-thick Chromium absorbing 
layer is shown in Figure 4-5(c) and (d) respectively. The amplitude and phase of the k-mask 
model is also shown for comparison. Finally, the two distinct aerial images resulting from 
the r-mask and k-mask models for imaging through a system with NA=0.7, R=4 and a=0.3 

are shown in Figure 4-5feV The 

w v ,* — .xww. liiw on. x oivi, cviucni oniy inrougii the 

r-mask model, is the pronounced difference of the two image lobes, with the lobe resulting 
from the 180° phase-shifted region having lower intensity than the lobe resulting from the 
unshifted 0° region. This kind of discrepancy between the two aerial images is usually 
attributed to "vertical" or "deep-trench" EM effects, that, as their cryptic names suggest, 
are impossible to model using a k-mask model. In the remaining of this Chapter, properly 
matched k-mask models, hereafter referred to as mk-mask models, capable of accurately 
including these effects in imaging are sought. 



Ln?v TC^f™ ""ft ^ "'""T ^ " 2 P,3ne *" t iS n0rmalIy incident on the P hotom ^ conse- 
where *S? t^TST ' nd,St ' nguishab,e - ™ e c " io " — here is the one noted in Chapter 2, 
where the TE has the E -field parallel to the y-axis and the TM parallel to the x-axis. 
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Figure 4-5. The intensity imbalance problem of a 0°/180° alt. PSM 

(a) Mask geometry, (b) amplitude of the near EM field solution from TEMPEST 
(X=193nm), (c) amplitude and (d) phase of the near EMF solution at the observation 
plane, located 50nm below the absorption layer, compared with the assumed k-mask 
model, (e) the resulting aerial images of the r-mask and k-mask models for an imag- 
ing system with NA=0.7, R=4, a=0.3. 
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4.2. Development of a new domain-decomposition and spectral-matchine 
methodology for alt. PSM with ID layouts 8 

Models for one-dimensional mask layouts are the subject of this Section and their 
treatment is essential for the development of models for the general case of two-dimen- 
sional layouts. 

tp£^l£*ta deCOmposition method (DDM > of an PSM into single- 

Consider the decomposition of an alt. PSM shown in Figure 4-6. This decomposi- 
tion is similar to the one shown in Figure 4-3(b) or [79], but it contains one extra term, 
namely that of the non-negligible transmission of the "opaque" Cr-based absorption layer,' 
whose importance was first noted by Yuan [115] and later incorporated in the model that 
Hotta used to study the trench design of alt. PSMs [46]. The pictorial representation of the 
decomposition method shown in Figure 4-6 implies that the scattered field below the mask 
can be adequately approximated by considering each mask opening by itself through dis- 
tinct rigorous mask simulations and properly accounting for the transmission of the absorp- 
tion layer that may have been counted more than once in the dark regions 1 . Such a 
decomposition method was seen to have its origin in the linearity of the Kirchhoff-Fresnel 
diffraction integral and predicts the scattered field with an accuracy of better than 0.5%, for 
a large number of practical situations, when compared with the true mask geometry with 
all openings present simultaneously. This decomposition was found accurate for shallow- 



§||^ well ~ jgF " tia ^^^^^| + |||'||||pl HI 





Figure 4-6. 0 l -order domain decomposition method (DDM) of an alt. PSM 

into single opening masks 



1 In the case of two openmgs (i.e. two single opening masks) the light below each dark region of the oriri 

£S5 " d0 H U ^ C °r*? Subtracti "S the uni *>™ transmission from the absorpt ive/ Us 8 
stra.ghtforward to venfy that in the case of an N-opening mask decomposed into N sinTle^em^ ™L 
the umform transmission of the absorption layer should be subtracted (N-l) times ^ 8 



72 



trench alt. PSM technologies, but it breaks down for deep-trench technologies because of 
cross-talk between adjacent features [46]. These limits will be explored in more detail in 
the following. The question of why one should bother decomposing a ID mask layout 
problem, that requires a fairly quick 2D rigorous simulation to determine the mask scatter- 
ing properties, into a set of ID masks, each requiring a 2D simulation, will be settled after 
considering the decomposition of masks bearing 2D layouts. 

For the 400nm (4X), 1:1 line/space pattern considered in Figure 4-5 the decompo- 
sition of the 0°/180° alt. PSM according to Figure 4-6 is shown in Figure 4-7. In Figure 4- 
7(a) the EM field of the original mask domain is decomposed into the first two single-open- 
ing domains shown in Figure 4-6. In Figure 4-7(b) the EM field across the observation 
plane is shown for both the original alt. PSM (both openings are present simultaneously) 
as well as the superposition of the two openings. The agreement between the two is almost 



(a) 




x - position (urn, 4X) real axis 



Figure 4-7. Example of the O tn -order DDM 

(a) Decomposition of the EM field of Figure 4-5(b) into the first two terms of 
Figure 4-6, (b) comparison of the near field after the decomposition shown in (a) 
with the near field of the exact alt. PSM geometry of Figure 4-5(a) at the observa- 
tion plane, (c) comparison of the "through-the-lens" (TTL) complex diffraction 
orders of the EM fields shown in (b). 
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perfect, but upon close inspection it is clear that the near fields differ slightly at the dark 
regions. This is because the third term of Figure 4-6 is also needed. Adding it gives an accu- 
racy better than 0.5% everywhere. Finally, Figure 4-7(c) depicts on the complex plane the 
three TTL (through-the-lens) diffraction orders (k k<,. k +1 ), for both the original alt. PSM 
and the decomposition shown in Figure 4-7(a). Note that the non-zero ko is responsible for 
the imbalance effect and that the small difference of the two methods in ko is fully attributed 
to the residual background transmission, which was not accounted for. 

4.2.2. Identifying the discrepancy between k-mask and r-inask models - The new mk- 
mask model 

Figure 4-8(a) shows the EM field snlntinn of a 

■ — — wpwuiug uuucf inpiane 

wave excitation. The amplitude of the rigorously calculated diffracted field at the observa- 
tion plane located 50nm below the absorbing layer is shown in Figure 4-8(b). Also shown 
in Figure 4-8(b) are the amplitude of the original k-mask model, which is just a unit-ampli- 
tude, 400nm-wide (4X) rect function properly centered, and the amplitude of the modified 
CMTF (mk-mask), which has resulted from a systematic procedure that will be explained 
in the next Section. For now, note that the modified (adjusted) mk-mask in this case is also 
a rect function with reduced width, slightly increased amplitude and slightly different phase 
than the k-mask. The result of this modified CMTF on the Fourier spectrum is depicted in 
Figure 4-8(c) and (d). The magnitude of all propagating diffraction orders vs. the wave- 
number k x is shown in Figure 4-8(c). Notice that the spectrum of the mk-mask matches the 
spectrum of the r-mask (rigorous CMTF) to better than 1% accuracy up to |k x |~0.012, 
which is well beyond the objective lens collection ability. Figure 4-8(d) illustrates nicely 
the effect of properly modifying the k-mask, where it is shown that the TTL diffraction 
orders can be relocated on the complex plane, so that they perfectly match (to within 0.5%) 
the orders resulting from the rigorously calculated r-mask. 
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Figure 4-8. Matching the k-mask of each individual opening to the r-mask 

(a) Electromagnetic field amplitude of a single 180° phase-shifted opening, (b) 
amplitude of the EM field at the observation plane shown in (a) (r-mask) com- 
pared with the original CMTF (k-mask) and the matched CMTF (mk-mask), 
.(c) magnitude of diffraction orders vs. k x and (d) the TTL complex diffraction 
orders on the complex plane for the three mask models shown in (b). 



A similar adjustment can be independently performed for a single 0° opening, 
which will result in a modified CMTF for that opening. Putting together the two mk-masks 
for the 0° and the 180° openings, as described in Figure 4-6, to form the mk-mask of the 
complete ID alt. PSM layout, results in a simulation of the aerial image that perfectly 
matches (to within <0.5%) the aerial image using the rigorous CMTF (r-mask). This is 
shown in Figure 4-9. 

The normalization used in spectra plots is such that a constant (DC) value of 1 
(clear field) would produce a Dirac 5-function at zero frequency with amplitude of 19.4 
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Figure 4-9. The mk-mask model for the example of Figure 4-5 

(a) Amplitude and (b) phase of the rigorous CMTF (r-mask) of the complete 
alt. PSM shown in Figure 4-5(a) and the matched CMTF (mk-mask). Note that 
the smaller opening corresponds to the 180° phase-well, (c) The TTL complex 
diffraction orders of the two CMTFs in (a), (b) are aligned and consequently 
the aerial images in (d) are identical. 



4.2.3. Systematic manipulation of the Fourier spectrum for the adjustment of the It- 
mask model 

In this Section we examine in depth how the attributes of the k-mask model can be 
adjusted to lead to a Fourier spectrum that closely matches the spectrum resulting from the 
r-mask. This is done through the example shown in Figure 4- 10(a), where a single 180° 
phase-well of CD=600nm (4X) with 50nm of isotropic glass underetch is considered. The 
mk-mask model (shown in Figure 4- 10(b)) is a properly centered rect function, but its 
width I^o, amplitude T 0 and phase <J> 0 are now varied until the best spectral match with the 
r-mask is achieved. The residual transmission of the background with amplitude T b and 
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phase <)> b are treated as constants, whose values are related to the thickness and the optical 
properties of the Chromium-based absorption layer. 
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Figure 4-10. The rect- function mk-mask model 



(a) Single 180° phase-well of CD=600nm (4X) with 50nm of isotropic underetch. 
The mask period is 4xCD=2.4um (b) The rect-CMTF model, where its width 
L x0 , amplitude T 0 and phase fa are varied until the best spectral match with the 
rigorous CMTF (r-mask) is achieved. T b and fa are the amplitude and phase of 
the residual background light. 



Figure 4-1 1(a) and Figure 4- 10(b) show the effect of reducing L xf> while keeping 
T 0 and <J> 0 constant. The spectrum of the r-mask is also shown for comparison. As expected 
from simple Fourier theory, reducing L x0 results in a broader amplitude spectrum 
(Figure 4-1 1(a)). The behavior of the seven TTL diffraction orders (up to +/- 3 for NA-0.8) 
on the complex plane as L x0 is varied is shown in Figure 4-1 1(b). 

Figure 4- 11(c) and Figure 4-1 1(d) depict the effect of reducing T 0 , while keeping 
L x0 and fa constant. Note that varying T 0 does not change the overall shape of the spectrum 
(the nulls in the amplitude spectrum of Figure 4-1 1(c) remain unchanged), but simply 
scales the amplitude level of the spectrum. 

Figure 4-1 1(e) and Figure 4-1 1(f) illustrate the effect of varying <J> 0 , while keeping 
L x0 and T 0 constant. The amplitude spectrum of Figure 4- 11(e) is completely unaffected 
by the changes in fa, but the phase of the diffraction orders is changing in liaison with fa y 
in such a way that increasing fa results in a counter-clockwise rotation of the spectrum, as 
shown in Figure 4-1 1(f) for the TTL diffraction orders. 
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Figure 4-11. Adjusting the mk-mask for good spectral matching 

(a), (c), (e) Magnitude of diffraction orders vs. k x of the rigorous CMTF (r-mask) 
and three rect-CMTFs (mk-masks), (b), (d), (f) the seven "through-the-lens" (TTL) 
diffraction orders (for NA-0.8). In (a), (b) T 0 =l, are kept constant and 

L x0 =600, 550, 500nm. In (c), (d) L x0 =600nm, 00=0° are kept constant and T^l, 
0.9, 0.8 and in (e), (f) L x0 =600nm, To=l are kept constant and $0=0°, 10°, 20°. 
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In summary, the adjustment of the three parameters of the matched Kirchhoff (mk- 
mask) model (L x0 , T 0 , %) proceeds as follows: First, L x0 is adjusted (usually reduced) 
from its originally designed mask value, CD=600nm (4X) in this case, until the overall 
shape of the k-mask amplitude spectrum matches that of the r-mask spectrum. A good mea- 
sure of the "degree of spectral matching" is to look at the location of the first nulls for the. 
two spectra. Next, T 0 (originally is 1) is raised or reduced accordingly, until the same TTL 
amplitude spectra are achieved. In practice, this model matched the amplitude spectra of 
all mask geometries that were considered well beyond |k x |=0.01, i.e. it is more than ade- 
quate even for ideal lenses of NA up to 1 . Finally, (|> 0 (originally set to the true phase of the 
EM field emerging from the phase-well with respect to an unetched clear background) is 
adjusted accordingly, so that the spectrum of the k-mask is rotated by the proper amount 
that will overlap it with the r-mask spectrum. 

The above steps have been automated through the use of custom routines written in 
the MATLAB environment. Although the aforementioned procedure for obtaining the mk- 
mask model is seemingly simple, there are nevertheless some points that require caution. 
Adjustment of L x0 , T 0 , (|>o should start from a judiciously chosen starting point. A "sharp" 
global optimum solution is usually not available, and, if the spectral matching is restricted 
only in the lens collection band, with L x0 , T 0 , § 0 allowed to take values in a large interval, 
the resulting solution for L x0 , T 0 , <t>o ma y be a set of values vastly different than the k-mask 
model. Although there is nothing fundamentally wrong with that, a solution that is "ele- 
gant" is preferred, in the sense that such a solution should have a "good degree of resem- 
blance" with the original k-mask. If such is the case, then the new mk-mask model makes 
more physical sense, since, in a way, it conveys information about the edge-diffraction 
effects that take place in an alt. PSM and that the k-mask model fails to capture. For exam- 
ple, as seen in Figure 4-8, the mk-mask for the 180° opening looks similar to the unadjusted 
k-mask, but it has a smaller L x0 , since owing to edge diffraction less light goes through the 
opening and the opening appears smaller than L x q. Its amplitude and phase are also 
adjusted so that they compensate for the average amplitude and phase distortions of the true 
r-mask for that opening. 
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4.2.4. Other possible functions for the mk-mask model 

Functions for the mk-mask other than the rect were also tried and evaluated. These 
included the Gaussian, the raised-cosine and one with local rect-functions near each edge 
to model the edge diffraction phenomena. They all yielded relatively successful spectral 
matching for bandwidths ranging from half to slightly more than that achieved using the 
above rect function mk-mask. Note that both the Gaussian and the raised-cosine can be con- 
structed to have good resemblance with the true scattered fields and subsequently offer per- 
formance equivalent to the rect-function. The model that attempted to replicate the edge- 
diffraction effects present to the mask geometry by adding local rect functions close to the 
edges was by 50% inferior (it only matched up to half of the bandwidth of the rect mk- 
mask). Although such a model might be intuitively satisfactory, and potentially versatile, 
since it suggests that regardless of the opening size the edge diffraction can be modeled in 
the same way, it was found to suffer from narrower bandwidth performance compared to 
the rect-CMTF, reminding us that edge-diffraction in the photomask is not localized near 
the edges, but rather extents over most of the feature opening. 

4.3. Towards a comprehensive quasi-rigorous method and the mk-mask 
for alt. PSM with arbitrary 2D layouts 

So far, an alternative methodology for the imaging simulation of alt. PSMs with 
arbitrary ID layouts has been developed. It consists of two basic steps: First, a decomposi- 
tion of the ID layout into multiple single-opening masks is performed, depicted in 
Figure 4-6, by virtue of the superposition principle of the Kirchhoff-Fresnel diffraction 
integral. This O^-order DDM was found accurate for all cases where cross-talk is negligi- 
ble. Next, through exhaustive rigorous simulation of each of these single-opening masks a 
matched-Kirchhoff-mask is found for each constituent and subsequently for the original 
layout. This compact mk-mask model, if used in the imaging simulation program, renders 
an image that is an excellent approximation to the r-mask model. The above methodology 
has no value, other than pedagogical, since it requires multiple rigorous simulations instead 
of just one with the complete ID layout in place. The same remarks are true in the extension 
of the 0*-order DDM in 2D layouts (Figure 4-3(a)). 
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However, the following conjecture can be used to revive the interest in the above 

methodology applied in 2D layouts: After the O th -order DDM has been performed in the 
2D layout, similar to Figure 4-3(a), but including the transmission of the absorbing layer 
as in Figure 4-6, no 3D rigorous simulation of the constituent masks is required! Instead, 
each constituent single-opening mask of the 2D layout is further decomposed as shown in 
Figure 4-12. According to this proposition, a quasi -rigorous result (qr-mask) of the scat- 
tered field can be obtained by executing two much faster 2D simulations and recombining 
the results as shown in Figure 4-12, where the last term is just a k-mask model that requires 
no rigorous simulation. There is no easy way to prove (or disprove) this conjecture. The 
large variety of cases that it works successfully, by giving accurate results is the closest 
thing to a proof. Remember also that our concern is an accurate calculation of the aerial 
images, and not near fields! Errors in near field calculations can be tolerated, provided that 
the aerial images are accurate. This is exactly what the k-mask model is, an inaccurate near 
field estimate that works well in aerial image calculations for certain situations. Note that 
the qr-mask model is expected to work well when the lateral dimensions of each feature 
are larger than X and cross-talk is negligible. 



(a) (b) (c) (d) 




1, 180° 1, 180° 1, 180° 1, 180° 



Figure 4-12. Illustration of the quasi -rigorous decomposition method 

Instead of the slow but completely rigorous 3D simulation of the original mask in 
(a), that would lead to the r-mask model, the scattered field below the mask is 
found by adding the r-masks shown in (b) and (c) and the k-mask shown in (d). 
Note that the r-masks in (b) and (c) require only 2D rigorous simulations that are 
much faster and the k-mask in (d) requires no simulation. The resulting approxi- 
mate mask model is dubbed the qr-mask (quasi-rigorous mask). 
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Next, two examples that demonstrate the validity and value of this method are pre- 
sented. The qr-mask model can either be used as shown in Figure 4-12, or each r-mask from 
the respective 2D rigorous simulation can be reduced to a compact mk-mask model, that 
matches the TTL-spectrum, as discussed in Section 4.2. 

4.3.1. Example 1: Dense line/space pattern 

Figure 4-13 depicts a 0°/180°, CD=100nm (IX), 1:1 dense line/space, simple layout 
of an alt. PSM with 50nm of isotropic underetch. The fully rigorous calculation of the aerial 
image of this mask under completely unpolarized illumination, with X=193nm, NA=0.7, 
o=0.3 and R=4 requires, according to the discussion in Chapter 2, the rigorous EM field 
solutions under both TE and TM plane wave excitations. Although a discretization of -25 
cells per X is sufficient for most practical cases [105], in this example the 3D mask domain 
was discretized with 40 cells per X, in order to achieve indisputable accuracy of the near 
EM fields of better than 0.5% and make the benchmarking comparisons with the quasi-rig- 
orous and mk-mask methods. Consequently the processing time and memory requirements 
for the 3D simulations rose dramatically. Each 3D TEMPEST simulation (TE and TM) ran 
for approximately 43hrs on a single 550MHz processor of the -100 CPU Millennium clus- 
ter at UC Berkeley [118] and used up to 1.57Gb of memory. Alternatively, the results can 
be obtained via the quasi-rigorous DDM and subsequent mk-mask models, as described in 
Sections 4.2 and 4.3. A total of eight 2D rigorous EM simulations are required that take a 
total of just under 12min to run and utilize negligible memory (~2.5Mb). 

The agreement of the mk-mask model with the full 3D rigorous EM field simulation 
(r-mask) is spectacular. The image intensity contours at 30% of the clear field virtually 
overlap for both TE and TM light and the same is true for the aerial images at the cut-lines 
near the center and near the edge of the lines, as shown in Figure 4-13(b) and Figure 4- 
13(c). For comparison, the aerial image contours and at the same cut-lines of the original 
k-mask (unmatched CMTF) without adjustment is also shown on these plots. 

Note that in this example there is virtually no cross-talk. The reasons for that will 
be obvious after the end of the next Chapter. 
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Figure 4-13. Examplel: Simple dense line/space layout 

(a) Mask layout, (b) aerial images under TE excitation and (c) aerial images 
under TM excitation. The 30% intensity contours are shown and the intensity 
profiles at the specified cut-lines. 
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4.3.2. Example 2: Semi-dense contact hole mask 

Figure 4-14 depicts a 90°/270°, CD=150nm (IX), 1:1.5 semi-dense contact hole alt. 
PSM with 50nm of isotropic underetch. The calculation of the aerial image of this mask 
under completely unpolarized illumination, with X=193nm, NA=0.7, a=0.3 and R=4 
requires again the rigorous EM field solutions under both TE and TM plane wave excita- 
tions. However, owing to the 90° rotational symmetry of this mask it suffices to run only 
one 3D simulation. Again, a discretization of the mask domain with 40 cells per X was used 
and TEMPEST ran for approximately 55hrs on a single 550MHz processor in the Millen- 
nium cluster and used 1 .82Gb of memory. A total of four 2D rigorous EM simulations are 
required by the qr-mask and mk-mask methods, that require a total of less than 8min to run 
and utilize approximately 3Mb of memory. The aerial images are compared with the one 
obtained using the conventional (unadjusted for edge-diffraction) k-mask model. Once 
more, the mk-mask model is able to capture accurately the edge-diffraction effects, since 
the image overlaps with the one obtained using the r-mask. The image with the k-mask 
model is seen to be vastly different. Cross-talk in this example is negligible. 
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CD=150nm(lX), 1:1.5 dense 
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Figure 4-14. Example 2: Simple semi-dense contact hole layout 

(a) Mask geometry of a 90°/270°, CD=150nm (IX), 1:1.5 semi-dense contact 
hole alt. PSM with 50nm of isotropic underetch. (b) Contours at 30% intensity of 
the aerial image under unpolarized light at A=193nm, imaged with NA=0.7, 
a=0.3 and R=4. Three images are compared, namely the one resulting from full 
3D rigorous simulation (r-mask), the one from a conventional scalar simulation 
with the k-mask model and the one from with the mk-mask model of the qr- 
DDM. (c) The same three images as in (b) at the cut-line that runs through the 
center of the contact holes. 
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43.3. Speed-up factor of the qr-DDM and mk-mask methodology 

In the examples presented, the speed-up factors of the qr-DDM and mk-mask meth- 
ods were seen to be -400. These speed-up factors would have been almost halved for a dis- 
cretization of -25 cells per X at the expense of reduced accuracy. The speed-up factor S 
resulting from the quasi-rigorous domain decomposition method (qr-DDM) combined 
with the mk-mask model is proportional to the size of the 2D mask layout, the discretiza- 
tion chosen for the rigorous EM field simulations and inversely proportional to the number 
of rectangles that constitute the layout of the dark-field alt. PSM (mask density): 

o a d 

Equation 4-4. a 2N 

where a is the length of the 2D layout (let length=width), d is the number of cells per X in 
the discretization of the mask domain and N is the number of rectangles that constitute the 
layout. The factor 2 enters because every rectangle requires two 2D simulations, according 
to the qr-DDM of Figure 4-12. 

Although speed-up factors of -400, while maintaining excellent accuracy, are quite 
impressive they are probably not sufficient for rapid aerial image simulations in an OPC 

environment, where the full-chip mask layout is broken up in coherence windows 1 and 
multiple small perturbations of the original mask layout within each coherence window 
should be evaluated in fractions of a second. To do this, another idea that is being used 
extensively in OPC software can be exploited, namely the idea of using look-up tables to 
further speed-up some calculations. To this end, a large variety of ID features that spans 
the whole interval of mask feature sizes that need to be addressed can be pre-calculated, 
converted to the compact mk-mask model and stored. Every time that a certain feature size 
is encountered in the layout, instead of using the k-mask model, the mk-mask model is 
retrieved from the look-up table for an accurate aerial image evaluation, since the edge- 
diffraction effects are properly modeled with the mk-mask. If the task of creating the look- 
up table is completed, then the notion of a speed-up factor is irrelevant, since the whole 



1. The notion of a coherence window can be found in [27]. It can be loosely defined as a mask fracturing, 
such that the image intensity at the center point of each fractured piece depends only on the layout within 
this piece. Proximity effects outside the coherence window can be safely neglected. 
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procedure is executed off-line and the speed is fixed - how fast the look-up table can be 
addressed. With mask feature sizes for every technology generation usually restricted 
within a relatively small set of values, the look-up tables are expected-to be fairly easy to 
generate, compact and consequently the look-up time short. 

A powerful alternative to the look-up table approach based on edge decomposition 
is presented in Chapter 6. 

4.3.4. Accuracy issues 

Two main concerns about the accuracy of the qr-DDM, and as a result the mk-mask 
model, can be identified: First, the corner effects are not modeled. Using the decomposition 
proposed in Figure 4-12, the solution for the scattered field below the mask is synthesized 
by two 2D simulations and as a result the electromagnetic interaction of light near the cor- 
ners is not accurately modeled. Specifically, the more strict boundary conditions that the 
field components need to satisfy near the corners are effectively substituted by a pair of 
more relaxed boundary conditions that hold near the edges of the absorbing layer in the 2D 
simulations. Second, the qr-DDM effectively does not allow for any cross-mixing of 
energy between the components of the scattered field to happen. If the illuminating light is 
polarized in the x-direction, then, according to the proposed method in Figure 4-12, all scat- 
tered energy will be polarized in the x-direction also. The same is true for the y-polariza- 
tion. However, cross-mixing of energy between field components is known to happen in 
real 3D masks and this cannot be modeled via the qr-DDM. 

It may seem puzzling by the two examples presented earlier how an aerial image 
was achieved that was hardly distinguishable by the one obtained with the fully rigorous r- 
mask model, albeit the accuracy concerns raised above. The reason for this can be under- 
stood if one considers the spatial extent of any corner effects and the spatial frequencies in 
the energy that is coupled into polarizations other than the excited one. The spatial extent 
of corner effects is typically much less than the wavelength. Therefore, it only affects the 
highest diffraction orders, that in a typical imaging system do not contribute at all to the 
image formation. A similar argument holds for the energy cross-coupled into polarizations 
other than the excited one, as seen in Figure 4-15. The mask from the example presented 
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Figure 4-15. Ignoring the perpendicular polarizations 

The example of Figure 4-14 is examined in more detail. The mask is illuminated with 
light polarized in the y-direction. Therefore most of the energy of the scattered near 
field is in the y-polarization, as seen from the graphs on the left. The non-negligible 
light in the x- and z- polarizations, due to cross-mixing, is seen to have a very high 
spatial frequency texture. This is exactly why, as seen at the plots on the right, the 
contribution of the x- and z- polarizations to the aerial image is negligible, at least by 
a factor of -300 smaller, (a) Amplitude of E x scattered field and contribution to the 
image, (b) and (c) the same for the E y and E z components respectively. 
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in Figure 4-14 is illuminated with an Ey polarization.However, some energy is seen to scat- 
ter in the and E z polarizations also. Although the local amplitude of the cross-coupled 
polarizations is seen to be up to 50% of the clear field, the spatial frequency of this light is 
so high that it does not contribute to the image formation. As seen in Figure 4-14, the 
second largest component in the image, which comes from the F^-polarized light, is smaller 
than the main Ey-polarized component by a factor of -300. Hence corner effects and cross- 
polarization effects can be safely neglected. 

4.3.5. The qr-DDM for arbitrary 2D layouts 

So far, the layouts considered were restricted to rectangularly shaped openings. In 
this Section the extension to arbitrary Manhattan-type polygons is developed. 

Consider the example layout shown in Figure 4-16. It consists of two unshifted rect- 
angular islands (with transmission of one and 0° phase) and a phase-shifted polygon of the 
Manhattan-type (with transmission of one and 180° phase). The natural extension of the 



1,0° 0 




Figure 4-16. qr-DDM of a layout into mutually disjoint rectangles (MDR) 
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domain decomposition method shown in Figure 4-3 is the one depicted in Figure 4-16, 
where the layout is decomposed into a set of mutually disjoint rectangles. Note that such a 
decomposition is not unique. Next, the quasi-rigorous conjecture of Figure 4-12 is applied 
on each rectangular layout component. This DDM is referred to as the qr-DDM with MDR 
(mutually disjoint rectangles). Alternatively, instead of a decomposition of the phase- 
shifted polygon into mutually disjoint rectangles, a decomposition into the set of maxi- 
mally overlapping rectangles (MOR) is performed as shown in Figure 4-17. The rectangu- 
lar intersections of the maximally overlapping rectangles are subtracted accordingly, so 
that no region of the original polygon is counted more than once. Note that this decompo- 
sition algorithm is particularly attractive not only because of the uniqueness of the set of 
rectangles that compose the polygon and have maximum overlapping area, but primarily 
because no artificial discontinuities are introduced. Observe in Figure 4-16 that, because 
of the "stitching" that occurs where the MDRs meet, the scattered field solution in those 
locations will be considerably different than the true solution of the complete polygon. This 
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Figure 4-17. qr-DDM of a layout into maximally overlapping rectangles (MOR) 
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is not the case in Figure 4-17, where the addition of the MORs and the subsequent subtrac- 
tion of their cross-sections is seen to remove all artificially introduced edges. Evidence of 
the validity of these arguments is given in Figure 4-18. The amplitude of the scattered Ey 
field as determined from a 3D rigorous solution of the complete layout under TE (Ey) exci- 
tation is shown in Figure 4- 18(a). In (b) and (c) the amplitude of the errors of the scattered 




(d) (e) (f) 




Figure 4-18. Comparison of MDR and MOR qr-DDM (1) 



(a) Amplitude of scattered E y field from 3D rigorous calculation of the whole 
mask (under E y excitation). The amplitude of the error of the scattered Ey field, if 
the mutually disjoint rectangle (MDR) or the maximally overlapping rectangle 
(MOR) qr-DDM is used, is shown in (b) and (c) respectively. Observe that the 
MOR method is superior to the MDR. (d) Aerial image of the mask shown in 
Figure 4-16 or Figure 4-17 using the r-mask model from (a). The amplitude of 
the error incurred in the aerial image simulation, if MDR or MOR qr-DDM is 
used, is shown in (e) and (f) respectively. Note that error in (e) approaches an 
unacceptable 10% of the clear Field whereas in (f) it peaks at 1%. Imaging system 
parameters used: X=193nm, NA=0.7, o=03 and R=4. 



E y field, if the MDR or MOR method of decomposition is used, is shown respectively. 
Observe the large error incurred near the "stitching" locations of the MDR method, which 
will lead to errors in the image simulation. On the other hand, the error of the MOR decom- 
position method is much smaller and concentrated locally near the "exterior" corners of the 
layout. Such corners are artificially constructed by the MOR method, as seen in Figure 4- 
17, and cannot therefore fully capture the true electromagnetic behavior of light. In spite 
of that, the spatially concentrated nature of the error implies that most of its energy is in 
the higher spatial frequencies that are discarded by the imaging system. Hence, they should 
not further contribute to errors in the imaging simulation. This is amply demonstrated in 
Figure 4-18 and Figure 4-19. The image of this layout if the r-mask model is used is shown 
in Figure 4-1 8(d) and the error of the MDR and MOR methods are shown in (e) and (f) 
respectively. The error level in the image is at an acceptable level of -1% for the MOR 
method but it locally approaches -10% for the MDR method. The images are further com- 
pared in Figure 4-19, where the 30% intensity contours are shown to virtually overlap for 
the r-mask and MOR method, but the contour of the MDR method differs. 




Figure 4-19. Comparison of MDR and MOR qr-DDM (2) 

The 30% contours of the aerial images resulting from the fully rigorous r-mask 
and the qr-DDM methods shown in Figure 4-16 (MDR) and Figure 4-17 (MOR) 
are compared. Note the excellent agreement between the r-mask and the qr-mask 
with MOR within the enlarged snapshot. The image resulting from the mutually 
disjoint rectangle qr-DDM is seen to differ substantially. 
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4.3.6. Off-axis illumination 

All simulation examples presented so far during the development of the domain 
decomposition methods involved normally incident light on the masks. This was because 
of the use of Hopkins' method for imaging, as described in Chapter 2. There is however 
nothing exceptional behind the use of normally incident light and the powerful decompo- 
sition techniques should be equally accurate for off-axis angles of incidence. In such an 
event, other imaging formulations, such as Abbe's, or a modified Hopkins' method that 
allows for the diffraction orders to depend on the angle of incidence, can incorporate sim- 
ilar decomposition methods that speed-up the imaging simulations without compromising 
the accuracy. 

In the example of Figure 4-20 a 0°/180°, 1:1 dense line/space alt. PSM with 
CD=400nm (4X) is illuminated by a plane wave that is traveling in air at -49° off-normal 
and owing to refraction is entering the glass of the mask (n glass (A,=193nm)= 1.563) at -29°. 
Note that in a typical imaging system the incident waves on the mask rarely exceed 10°- 
15°. The r-mask simulation across the observation plane is compared with the qr-DDM of 
the mask into the two openings. The two curves are identical! Observe that the scattered 
fields are highly asymmetrical and mk-mask models should be adjusted accordingly. The 




Figure 4-20. The qr-DDM under highly oblique illumination 

The qr-DDM is seen to be accurate even for off-axis plane wave illumination at 
very large angles of incidence. A 0°/180°, 1:1 dense line/space alt. PSM with 
CD=400nm (4X) is illuminated by a plane wave incident on the mask at 0ai r ~49° 
(6 gIass ~29°). The amplitude of the scattered field along the observation plane is 
compared for the r-mask and the qr-DDM. 
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simple rect function mk-mask model will not perform adequately for such high angles of 
incidence. 

4.4. Conclusions 

Rigorous electromagnetic simulation of the role of the mask in the imaging process 
has so far been limited to 2D -problems and relatively small periodic 3D problems. 
Although this capability is enough to tackle a variety of advanced photomask design issues, 
it is unquestionably restrictive. For example, rigorous simulation of a lOjim by lOjum area 
on an alt. PSM that contains random (non-periodic) structures, at A,=193nm, would require 
on the order of ~3.6Gb of memory (discretization of the ~14,000A. 3 domain with 20 cells 
per X) and would run for many hours on a powerful CPU. Although the scalar theory can 
handle large mask domains, it has been discarded for advanced photomask technology sim- 
ulations as being inaccurate. 

In this Chapter, domain decomposition methods were developed that allow the 
division of a larger mask diffraction problem into a set of constituent parts, followed by the 
application of simple models for each constituent part that match the diffraction spectrum 
and subsequent synthesis of the original problem using these simplified models of the con- 
stituent parts. These methods originate from the linearity of the Kirchhoff-Fresnel diffrac- 
tion integral. A ID layout is decomposed into a set of single-opening masks and 
subsequently each single-opening mask is treated separately. Furthermore, a compact rect- 
function model for the diffracted near field of each opening is extracted, such that it 
matches the through-the-lens (TTL) spectrum of the rigorously calculated near field (r- 
mask). Application of the method in 2D layouts is based on Figure 4-12 and the conjecture 
that only 2D rigorous electromagnetic simulations will be used to synthesize diffraction 
from 3D structures. The spectral matching method is again used to simplify the r-mask 
models of each constituent 2D rigorous electromagnetic simulation. The excellent accu- 
racy of the method applied in 2D layouts (3D diffracting structures) was attributed to the 
fact that errors arising from neglecting corner effects and cross-coupling of energy into 
light polarizations other that the excited one, are mapped at the extremities of the spectrum 
of propagating plane waves and do not contribute to the image formation. The decomposi- 
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tion of arbitrary 2D (Manhattan-type) layouts was shown to be accurate if it follows the rule 
of "maximally overlapping rectangles", whereas the decomposition into "mutually disjoint 
rectangles" is less accurate. Finally, it was shown that the decomposition method is valid 
even when off-axis illumination is used. 

The excellent accuracy (better than 99% compared to the r-mask models) in the 
examples presented in this Chapter was accompanied by speed-up factors exceeding 200. 
It was stressed though that an off-line pre-calculation of the diffraction from various 2D 
geometries and creation of a look-up table of mk-mask models can further speed-up the cal- 
culations. The benefits of such accuracy-speed combinations can be reaped for example by 
model-based OPC software tools, that so far rely on the inaccurate scalar theory for aerial 
image calculations. Conceivably, accuracy on par with a rigorous 3D simulation becomes 
feasible for large - even full mask - areas. 

The expansion of the domain decomposition methods of this Chapter to include the 
complex effects of energy cross-coupling between neighboring mask apertures is treated in 
the next Chapter and a different (more powerful) approach of domain decomposition based 
on edges is developed in Chapter 6. 
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5 Deep Phase- Well and Cross-Talk 
Effects in Alternating Phase-Shifting 
Masks 



Simulation can give insight into the complex electromagnetic effects of edges 
and cross-talk between edges. In this Chapter the domain decomposition technique is 
expanded in order to isolate and quantify these phenomena. 

In order to avoid confusion in the subsequent treatment the following nomenclature 
is adopted: Referring to the schematic diagram of an alt PSM in Figure 5-1, the glass- 
etched regions are the two phase-wells, namely the "deep" and the "shallow" phase- well, 
which always have an etched-depth difference that results in 180° phase-shifting of the EM 
field emanating from each opening. Each phase-well for a rectangularly shaped mask 
opening typically resembles a parallelepiped in three dimensions and as such it comprises 
of 6 faces. At the bottom face of each phase-well 4 edges and 4 corners form. Depending 
on the details of the actual mask fabrication process the edges and comers may be rounded, 
whereas the side-faces may be sloped (non-vertical). The two dimensional schematic of an 

edges 




side-faces (sidewalls) 
Figure 5-1. Schematic diagram of an alt. PSM 
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alt. PSM in Figure 5-1 shows only a cut-plane, hence the corners are not shown, the edges 
are depicted as corners and the faces are depicted as edges. 

5.1. Towards a comprehensive model for cross-talk 

Qualitative arguments that explain the cross-coupling of energy between neighbor- 
ing openings in alternating phase-shifting masks have been proposed by other researchers 
in the past. Socha [84] suggested that the presence of induced currents on the absorbing 
layers of a mask (with non-zero imaginary part of the refractive index) results in a "shorting 
out" of the transmitted fields near the boundaries of each opening and consequently 
decreases the amount of light passing through the openings. Wong [106] observed the ten- 
dency of the electromagnetic field to concentrate in the higher refractive index glass, as 
opposed to the air regions, in the vicinity of an air/glass interface, such as the faces of 
phase- wells and attributed the differences in the transmission of light between the two adja- 
cent openings (imbalance effect) in this behavior. Levenson et al [50] recently proposed a 
novel design for alternating PSM. Among other things, it successfully addresses the issue 
of cross-coupling of energy between adjacent phase-wells. In the following Section, a new 
model for the representation and calculation of cross-talk is introduced. 

5.2. Scattering off of a 90° air/glass discontinuity 

Figure 5-2 depicts the rigorous EM field solution for the instantaneous field being 
established, when a 90° air/glass discontinuity is illuminated from the top by a linearly 
polarized, monochromatic and normally incident plane wave. The bottom-left part of the 
plot is the air region with refractive index n=l and the remaining part, indicated by the 
dotted line, is the glass-filled region with refractive index n=1.563, at X=193nm. Note the 
inherent with the simulation program periodicity along the x-axis, which unavoidably leads 
to simulating a periodic array of 90° air/glass discontinuities. However, this periodicity 
does not pose a significant limitation, since the dimensions chosen are sufficiently large 
compared to X to mimic an isolated discontinuity. Observe from Figure 5-2(a) that the 
reduced phase velocity of light in glass causes the wave in the glass region to lag compared 
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to the wave in the air region, and also the fact that there exists a finite "dark" region to the 
left and to the right of the air/glass interface, in agreement with Wong's findings [106]. 

(a) (b) 




Figure 5-2. Scattering off of a 90° air/glass discontinuity 

(a) Instantaneous EM field being established when a 90° air/glass discontinuity is illu- 
minated from the top by a normally incident TE plane wave at k=193nm, (b) the 
instantaneous scattered field (at the same instant as (a)) because of the 90° corner. 



In order to gain more insight on the effect of the corner itself, the instantaneous 
field (at the same instant!) of a glass/air interface with no discontinuity along the x-axis 
was subtracted from the left-half of the original solution and the instantaneous field of a 
uniform plane wave travelling in glass was subtracted from the right-half of the original 
solution. The result of this process is the strange looking plot of the scattered field from the 
corner, depicted in Figure 5-2(b). There are two interesting observations to be made on this 
plot: First, a 90° corner at an air/glass interface results in a scattered wave that is cylindrical 
in nature with a non-uniform radiation pattern, where most of the energy is directed in the 
forward direction (direction of the incident field). Second, inside the glass at the right side 
of the air/glass interface, but below the location that the instantaneous incident field has 
had time to propagate there exists a region of non-zero scattered field. Clearly, the source 
of this field cannot be the corner discontinuity itself, since the scattered field from the 
corner has not had time to directly propagate in glass to that point. In order to satisfy the 
boundary conditions of Maxwell's equations everywhere along the vertical interface, a 
non-zero field should exist in the glass side, although there is no time for either the incident 
or the scattered field to reach that location. In a way, the "speedier" field on the air side is 
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depositing elementary Huygens' spherical sources [37] that radiate into the glass side. 
Since the velocity by which these Huy gens' sources are deposited is faster than the phase 
velocity of the disturbance that they cause inside the glass region, the end result is a "shock 

wavefront" of Cherenkov radiation 1 [47]. The angle by which this plane wave is tilted is 
theoretically given by: 



which is in excellent agreement with that found from our simulation. This leaking field is 
responsible for energy being "sucked" into the higher refractive index material and is also 



the ripples present in the plot of Figure 5-2(b). Yet another way to explain this field behav- 
ior has its origin in the phenomenon of total internal reflection [7]. Total internal reflection 
results when an EM field that is traveling inside a higher refractive index material (in this 
case glass with n= 1.563) impinges at an oblique enough angle on the interface with a lower 
refractive index material (in this case air with n=l). The boundary conditions along the 
interface impose this bizarre behavior, where no energy is reflected backwards in the opti- 
cally denser (higher refractive index) medium and all the energy is "trapped" in the form 
of surface waves travelling along the material interface [7]. In the case of the scattering off 
of the 90° air/glass discontinuity this phenomenon is observed in reverse, whereby the nor- 
mally incident plane wave along the air/glass interface oriented in the z-axis plays the role 
of the surface wave, which causes the field on the glass, side to take off at an angle equal to 
the critical angle. Sure enough the critical angle in this case is 39.78°. 

Next, the model for cross-talk between phase-wells is developed. 

5.3. Modeling the cross-talk between phase-wells 

The wave behavior at the air/glass discontinuities (faces, edges and corners) of the 
etched phase-wells of alt. PSMs suggests the model for cross-talk that is depicted in 



1. Strictly speaking, Cherenkov radiation involves relativistic motion of elementary particles. Nevertheless, 
the physical phenomenon that is observed here follows the same principle. 
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Figure 5-3. Figure 5-3(a) shows the geometry of a deep-trench alt. PSM with a dense line/ 
space pattern. The incident illumination causes all edges (corners in the 2D schematic dia- 
gram) present to scatter a non-uniform in angle, cylindrical radiation pattern, in a way sim- 
ilar to the result of the previous paragraph. Using a geometric ray tracing approach 1 and 
considering the reflection and transmission coefficients at all interfaces present (glass/air, 
air/glass, glass/absorption layer) one can identify "favorable angles" in the radiation pat- 
tern of the edges that do not contribute to cross-talk, shown in Figure 5-3(b) and also "unfa- 
vorable angles" that contribute to cross-talk, such as those shown in Figure 5-3(c) and (d). 
Note that, according to this model and the development that follows, the edges at the 



(a) (b) 




shallow phase-well 



(c) (d) 
unfavorable angles (cross-talk) 




Figure 5-3. A comprehensive model for cross-talk 

(a) Geometry of a deep-trench alt. PSM with a dense line/space pattern. 
All edges at the bottom of the phase- wells scatter radiation non-uniformly 
in all angles, (b) Cones of scattered radiation that do not contribute to 
cross-talk, because scattered light is reflected backwards or at highly 
oblique angles, (c) "Direct" cone of scattered light - cross-coupled light 
exits at a wide range of spatial frequencies, (d) "Reflected" cone of scat- 
tered light - cross-coupled light exits mostly at near-normal incidence. 



1- The dimensions of the problem are on the order of a wavelength. The geometric optics approach is 
not accurate in this regime, but it is nevertheless invoked to offer a qualitative explanation of the argu- 
ment. 
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bottom of the phase wells are identified as the primary contributors to cross-talk and nei- 
ther the faces, nor the corners play a significant role in this effect. This proposition is quite 
different from the prevailing ideas that the faces of the phase-wells play the most significant 
role in cross-talk and that the corners and the edges are also equally important. This model 
will be validated from its successful application in predicting and calculating cross-talk in 
a wide range of situations. 

5.4. 1 st order domain decomposition method of alt. PSM 

The zero-order domain decomposition of alt. PSMs presented in the previous Chap- 
ter is revised in order to include the cross-talk between adjacent phase-wells. This (first 
order) decomposition method is depicted in Figure 5-4 for iD layouts and contains two 
extra terms, namely CT 12 and CT 21 , that represent the amount of energy coupled from the 
shallow phase-well to the deep phase-well and the deep to the shallow respectively. 



Figure 5-4. l st -order domain decomposition method (DDM) of an alt. PSM 

The calculation of the cross-terms CT 12 and CT 21 can be performed using the sim- 
ulation experiment that is pictorially depicted in Figure 5-5. In order to determine the por- 
tion of scattered light from the edges of the shallow phase-well that goes through the 
opening of the deep phase-well, one proceeds as depicted in Figure 5-5(a): First, a simula- 
tion is run, in which the opening of the shallow phase-well is shut with the absorption layer 
but the phase-well itself is left intact. Then, from the resulting transmitted EM field across 
the observation plane of this geometry, the EM field solution coming from a simulation in 
which only the single opening of the deep phase-well is present is subtracted. Finally, the 
two remaining simulations shown in Figure 5-5(a) have to be run and the transmitted EM 
fields need to be subtracted/added, owing to the finite residual transmission of the absorp- 
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tion layer. Similar simulation steps are followed in order to calculate CT 2l . The cross- 
terms CT 12 and CT 2 i can be viewed as perturbation terms necessary to restore the accuracy 
of the zero-order DDM, since for openings that are spaced far apart these terms are seen to 
converge to zero. The total number of two dimensional simulations needed in the first- 
order decomposition is four instead of two for the zero-order decomposition (because of 
cancellation of terms), so the expected increased accuracy of this decomposition costs a 
doubling of required simulation time. 

Two specifically chosen examples of the first order decomposition are presented in 
a parallel fashion next and they demonstrate both the accuracy of the method and also 
resolve the misconception that cross-talk in alt. PSMs always leads to poorer aerial images 
with larger imbalance between the adjacent phase-shifted openings. 

A 1:1 dense line/space geometry with a mask CD of 400nm is patterned on a deep- 
trench 180°/360° alternating PSM in two ways: The first mask has no undercut and per- 
fectly vertical sidewalls (faces) and the second mask has 50nm undercut and sloped faces 



(a) 



CT 12 = 





CT 21 = 



(b) 





Figure 5-5. Calculation of the cross-terms: 

(a) CT 12 represents the portion of scattered light from the corners of the shallow 
phase- well that goes through the opening of the deep phase- well and (b) CT 2 \ repre- 
sents the portion of scattered light from the corners of the deep phase- well that goes 
through the opening of the shallow phase-well 
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as shown in Figure 5-6. Both masks are illuminated by a normally incident (travelling along 
the z-axis) plane wave that is TE linearly polarized (E field // y-axis) with A=193nm. The 
cross-talk terms CT 12 and CT 21 for each mask are calculated according to Figure 5-5 and 
their amplitudes are plotted in Figure 5-6(a) and (b). Next, the near diffracted fields from 
each mask along the observation plane are calculated in three ways: First, by simulating the 
complete mask structure with both openings present (r-mask), then by decomposing the 
mask using the zero-order domain decomposition method (Figure 4-6) and finally using the 
first-order domain decomposition that models cross-talk. The amplitude of these near dif- 
fracted fields is plotted in Figure 5-6(c) and (d) and the error of both decomposition meth- 
ods is plotted in Figure 5-7. Clearly, in this situation, ignoring the cross-coupled energy in 
the domain decomposition would be erroneous. The agreement of the first-order decompo- 
sition with the solution from the exact geometry is spectacular. The amplitude of the error 
is everywhere less than 2% of the clear-field value (which is 1), instead of as much as 
-30%! The normalized mean square error 1 (NMSE) for vertical sidewalls drops to 3.9E-4 
from 1.5E-2 when cross-talk is included. For the sloped sidewalls alt. PSM it also drops 
from 3.2E-3 to 2.4E-4. A key observation at the plots of Figure 5-6(c) and (d) is the fact 
that cross-talk in the case of the vertical sidewalls PSM is pushing the peak of the electric 
field emanating from the shallow phase- well to a lower value than had cross-talk been neg- 
ligible, whereas the opposite happens in the case of the PSM with sloped sidewalls, where 
the peak value of the shallow phase-well increases. The result of this is that in the first case 
cross-talk is helping balance the aerial image, whereas in the second case it negatively 
affects the balance, as shown in Figure 5-8 for an imaging system with NA=0.7, a=0.3 and 
R=4 at best focus. From these graphs it is also apparent that neglecting cross-talk leads to 
inaccurate calculation of the aerial images. The peak intensity % error of each lobe and the 
edge position error for a simple threshold resist model at 30% intensity are indicated on the 
plots. Note the excellent agreement in the aerial images after the inclusion of the cross-talk 
terms, where the images virtually overlap. 



1. The normalized mean square error is defined by: NMSE = j\E r - E^dS/^\E^dS , where E r is 
the "true" field resulting from a rigorous simulation of the complete mask geometry and E s is the syn- 
thesized field via the l st -order-DDM. 
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Figure 5-6. Examples of l st -order domain decomposition method: 1:1 dense line/ 
space pattern, CD=400nm (4X), 180°/360°, A=193nm, two different sidewall 

profiles are considered 

Normalized amplitude of CT l2 and CT 2 i for vertical (a) and sloped (b) sidewall 
profile. Comparison of the three simulation approaches: Complete alt. PSM sim- 
ulation (no decomposition, r-mask), 0 th -order and l st -order decomposition for 
vertical (c) and sloped (d) sidewall profile. 
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Figure 5-7. Errors of the 0 th - and l st -order decompositions for the examples of 

Figure 5-6 

Normalized amplitude of errors for vertical (a) and sloped (b) sidewall profiles. 
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Figure 5-8. Aerial images for the alt. PSMs of Figure 5-6 

System parameters: NA=0.7, (T=0.3, R=4, best focus. Normalized image inten- 
sity for vertical (a) and sloped (b) sidewall profiles. Again, the three simulation 
approaches are compared: Complete alt. PSM simulation (no decomposition, r- 
mask), O^-order and l st -order decomposition. 



5.5. Cross-talk dependence on the geometrical characteristics of the alt. 
PSM 

Using the above model, a systematic evaluation of the dependence of cross-talk on 
the geometrical details of the alt. PSM is performed. Although only the results for the terms 
CT 12 and CT 21 are presented, in each case the complete alt. PSM geometry was simulated 
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separately and the true transmitted field (r-mask) was compared with the field obtained 
through the first-order decomposition. The local error incurred was never larger than 2% 
of the clear field value. 

Figure 5-9 shows the cross-talk of four possible alt. PSM technologies, namely 0°/ 
180°, 90°/270°, 180°/360° and 270°/450°, with vertical sidewalls under TE illumination 1 . 
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Figure 5-9. Dependence of cross-coupling on phase-well depth 

(a) Geometry, (b), (c) normalized amplitude of cross-talk terms 
CT 21 and CT 12 for 0°/180°, 90°/270°, 180/°360° and 270°/450° alt. 
PSMs, (d), (e) normalized magnitude of CT 2 \ and CTj 2 . 
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Both the near field correction (cross-talk) terms CT J2 and CT 21 are shown along with their 
respective spectra in each case. It is worth to note that the deeper phase-well is contributing 
more cross-talk to the shallower (i.e. CT 21 >CT 12 ) in every technology and that the cross- 
talk increases for deeper phase-well technologies. Both these results are intuitively satis- 
factory, since according to the physical understanding of cross-talk based on Figure 5-3 the 
edges at the bottom of deeper phase-wells are located at more advantageous positions to 
throw their scattered light into the adjacent shallower phase-wells. Observe also from the 
spectra plots that most of the cross-talk energy is in specular directions and does not go 
through the optical system. This can be further quantified by considering the normalized 
total cross-talk energy and normalized through-the-lens (TTL) cross-talk energy, defined 
by Equation 5-1 and Equation 5-2: 



Equation 5-1. 5 = 

CT- total 



Equation 5-2. 5 _ 1*1 g k RW 

CT- TTL 



l*U*o 

J CT'ik^dk^ 



j ?{k x ,k y )dk x dk y 



The bar charts of Figure 5-10 show the percentage of total and through-the-lens 
(TTL) cross-talk energy for the 0°/180°, 90°/270° 180°/36O° and 270°/450° technologies. 
The lower portion of each bar (blue color) corresponds to the cross-talk energy from the 
shallow to the deep phase-well and the upper portion (red-brown color) to the cross-talk 
energy from the deep to the shallow phase-well. Observe that indeed the through-the-lens 
portion of cross-talk energy in approximately one order of magnitude smaller than the total 
cross-talk energy, which indicates that most of the cross-talk energy is directed in highly 
oblique angles and is lost in the sidewalls of the projector column. 



1. TM excitation causes similar cross-talk profiles with slightly reduced peak values. 
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Cross-talk vs. phase-well depth 
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Figure 5-10. Cross-talk energy dependence on phase-well depth 

Next, it is important to determine the dependence of cross-talk on pitch (period), 
since the present model points out that cross-talk should be insignificant when the lateral 
separation of adjacent phase-wells is relatively large, in which case the simpler zero-order 
decomposition method will suffice for the simulation of mask layouts. Figure 5-1 1 depicts 
the total and through-the-lens normalized cross-talk energy for a 180°/36O° alt. PSM with 
a fixed line of CD=100nm (IX) and variable pitch. Three values of pitch are shown, 
namely 200hm, 250nm and 400nm, which correspond to 1:1 dense line/space, 1:1.5 semi- 
dense line/space and 1:3 sparse line/space patterns. The key result is that indeed the TTL 
cross-talk drops quickly with increasing pitch. This is however an artifact of the normal- 
ization. The larger pitches are a result of larger openings (spaces). On one hand the abso- 
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Figure 5-11. Cross-talk energy dependence on pitch (mask period) - fixed line-width 
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lute amount of cross-talk does not drop, since the relative position of edges of adjacent 
phase-wells remains constant - although the period increases, but on the other hand the 
larger openings increase the total transmitted energy making the normalized cross-talk 
smaller, since the denominator of Equation 5-2 is larger. Therefore one must be careful 
when to discard the first-order domain decomposition method in favor of the zero-order 
when dealing with large pitches that are a result of large openings. 

The situation is different when the larger pitch is a result of increasing the line- 
width while keeping the space-widths (openings) constant. In that case the relative separa- 
tion of edges of adjacent phase-wells increases, the unfavorable light paths that cause cross- 
talk according to the model of Figure 5-3 are less and ™ * ~ c »it tu~ ,Ko n i„. . <- 

_ — — tuv uujwiuit OllIUUIIL or 

cross-talk energy quickly drops. This is evident by the bar charts of Figure 5-12 that show 
the total and through-the-lens normalized cross-talk energy for a 180°/360° alt. PSM with 
fixed space-widths of CD=100nm (IX) and variable pitch. Again, three values of pitch are 
shown, namely 200nm, 250nm and400nm, which correspond to 1:1 dense line/space, 1.5:1 
semi-dense line/space and 3:1 sparse line/space patterns. 
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Figure 5-12. Cross-talk energy dependence on pitch (mask period) - fixed space-width 

The cross-talk dependence on the amount of undercut (underetch) is shown in 
Figure 5-13, which depicts the normalized cross-talk energy for Onm, 50nm and lOOnm of 
undercut, for a lOOnm (IX) 1:1 dense line/space, 180°/360° alt. PSM with vertical side- 
walls. Again, in order to explain the drop of cross-coupled energy with increasing undercut, 
the model of Figure 5-3 can be recalled, according to which the increased undercut pulls 
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the faces of the phase- wells behind the opening in the absorber layer and consequently can- 
cels some of the light-paths that contribute to cross-talk. 
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Figure 5-13. Cross-talk energy dependence on undercut 

Figure 5-14 depicts the cross-talk dependence on the sidewall angle of the faces of 
phase-wells. The lOOnm (IX) 1:1 dense line/space, 180°/360° alt. PSM example is chosen 
again to demonstrate the effect. The three bars on each chart correspond to 90° vertical 
sidewalls, -80° sidewalls (labeled "tilted") and -70° sidewalls (labeled "more tilted"). The 
behavior is now somewhat erratic, with the cross-talk energy appearing unchanged for 
slopes up to -80° and then spiking up for larger sidewall angles. Since the integrity of the 

sidewalls in the alt. PSM fabrication process is typically good and very close to 90°, cross- 
talk changes are insignificant for such small sidewall slopes. 
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Figure 5-14. Cross-talk energy dependence phase-well profile (sidewall angle) 
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Finally, owing to the etching process in fabricating an alt. PSM sometimes the 
edges at the bottom of the phase-wells are not perfectly sharp, but they are rounded. The 
effect of the edge sharpness on cross-talk is shown in Figure 5-15, for the lOOnm (IX), 1:1 
dense line/space, 180°/360° alt PSM. The bar labeled "rounded" shows the cross-talk when 
the edges have a 50nm curvature. Cross-talk clearly increases with increasing edge curva- 
ture. This is almost counterintuitive, since one might be inclined to think that a rounded 
edge scatters less light. However, the also rounded edge from the neighboring shallow 
phase-well has now a better "acceptance ability" of the light scattered from the edges of the 
deep phase-wells and consequently more cross-coupling occurs. 



Cross-talk vs. rounded edges of phase-wells 
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Figure 5-15. Cross-talk energy dependence phase-well profile (rounded edges) 



5.6. Ideas for cross-talk elimination 

Based on the physical understanding developed for the origin of cross-coupling 
between adjacent phase-wells, two ideas for elimination of cross-talk are examined. Both 
these ideas also provide an indirect validation of the cross-talk model presented in 
Section 5.4. Note that eliminating cross-talk does not automatically solve all the intrinsic 
difficulties in controlling the details of the diffracted fields from the two adjacent phase- 
shifted openings. For example, cross-talk was shown in Figure 5-6 and Figure 5-8 to either 
suppress or increase imbalance. Eliminating cross-talk definitely makes the analysis easier, 
since there is one less parameter to worry about, but at the same time a potential degree of 
freedom for engineering the diffracted fields is lost. 
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Covering the faces (sidewalls) of phase-wells with the absorbing layer blocks both 
the direct light paths of Figure 5-3(c) and the reflected light paths (reflected from the back 
side of the absorbing layer) of Figure 5-3(d) and should consequently completely eliminate 
cross-talk. This is demonstrated in Figure 5-16 that depicts the amplitude of the cross-talk 
terms for a lOOnm (IX), 1 : 1 dense line/space, 180°/360°, vertical faces alt. PSM before and 
after covering the faces with a 50nm Cr-layer. After covering the sidewalls cross-talk drops 
to insignificant levels of less than 2%. This idea was implemented by Levenson et al. [50] 
as a means of controlling the imbalance and improving the performance of alt. PSMs for a 
wide range of pitches (periods). 



k=193nm 
E y (TE) 



welll 



80nm 




T 



r-t 

400nm (4X) 



1:1 line/space 
CD=100nm (IX) 



0.3 

0) 
N 

S 0.2 

§ 



0.1 



< 







exposed sidewalls 

— Cr-covered sidewalls 




1 — i I < 





0.3 

<D 
N 
-» 

£ 0.2 
O 

c 



exposed sidewalls 

— Cr-covered sidewalls 



0.4 0.8 1.2 

x - position (jwm,4X) 



1.6 




0.4 0.8 1.2 

x - position (fxm,4X) 



1.6 



Figure 5-16. Cross-talk elimination through the use of "blinders" 

The sidewalls are covered with a 50nm thick layer of Cr and the cross-cou- 
pling is eliminated. 



One of the advantages of simulation is that it can evaluate situations where exper- 
iments are too difficult to control. To this end, instead of fully covering the faces of both 
phase-wells, the shallow phase-well was left intact (no absorbing layer deposited on its 
faces) and the deep phase-well was covered half-way with absorbing layer as shown in the 
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schematic diagram of Figure 5-17. According to the physical understanding and model of 
cross-coupling developed in Section 5.3, the scattered light from the edges at the bottom of 
the shallow well cannot go through the opening of the deep phase-well and consequently 
this component of cross-talk is eliminated as shown in Figure 5- 17(b). On the other hand, 
the scattered light from the edges at the bottom of the deep phase-well can still go through 
the opening of the shallow phase- well unobstructed and therefore this component of cross- 
talk remains almost unchanged as shown in Figure 5-17(a). Note that the small increase of 
this cross-talk term is attributed to the extra scattered light caused by the discontinuity of 
the half-covered sidewalls of the deep phase-well with a 50nm thick layer of Cr. This sim- 
ulation successfully demonstrates the validity of the physical model in place for cross-talk. 
Had the primary source of cross-talk been the whole face (sidewall) of phase-wells and not 
the edges at the bottom, the amount of cross-talk from the deep to the shallow phase-well 
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Figure 5-17. Validation of cross-talk model 



The faces of the deep phase-wells are half-covered with 50nm of Cr-Iayer. This 
reduces dramatically the cross-talk from the shallow to the deep phase-well, whereas 
the cross-talk from the deep to the shallow remains almost unchanged. 
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would have been almost halved, when half of each face is covered, and would have not 
remain almost unchanged as this simulation demonstrated. 

The second idea for cross-talk elimination is based on blocking the reflected light 
paths (reflected from the back side of the absorbing layer) of Figure 5-3(d). In order to do 
that, a 50nm thick layer of M0O3 (with a refractive index of n=1.6+j0.69 at A=193nm) was 
deposited before the 80nm thick absorbing Cr-layer that acts as an anti -reflective coating 
(ARC). This anti-reflective behavior is shown in Figure 5-18, where the addition of the 
50nm M0O3 layer is seen to significantly reduce the amplitude of the standing wave being 
established behind the 80nm thick Cr-layer to levels that equalize it with the standing wave 
amplitude for a glass/air interface. The amplitude of the cross-talk terms for a lOOnm (IX), 
1:1 dense line/space, 180°/360°, vertical faces alt. PSM without and with the 50nm M0O3 
ARC is shown in Figure 5-19. The cross-talk levels with the ARC have been suppressed 
by almost 50%. This is indicating that both the direct light paths of Figure 5-3(c) and the 
reflected light paths (reflected from the back side of the absorbing layer) of Figure 5-3(d) 
contribute approximately equally to the cross-coupling effect in alt. PSMs. 
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Figure 5-18. A 50nm thick layer of M0O3 acts as anti-reflective coating 
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Figure 5-19. Cross-talk elimination through the use of ARC 
The cross-talk levels reduce by -50% when a 50nm M0O3 layer is used as ARC 
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Figure 5-20. Aerial images for the alt. PSMs of Figure 5-16 and Figure 5-19 



System parameters: NA=0.7, a=0.3, R=4, best focus. Normalized image inten- 
sity for Cr-covered sidewalls (a) and use of 50nm Mo0 3 ARC (b). The three sim- 
ulation approaches are compared: Complete alt. PSM simulation (r-mask, no 
decomposition), 0 th -order (no cross- talk terms) and l st -order domain decomposi- 
tion (including cross-talk). 
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The plots of Figure 5-20 depict the intensity of the aerial image for the alt. PSM 
geometries with the Cr-covered sidewalls and the 50nm Mo0 3 ARC. Inclusion of the 
cross-talk terms in the calculation is shown to be redundant, since the reduced level of 
cross-coupling in both cases makes the images obtained from the simulation of the exact 
geometry (r-mask) and from both the zero- (no cross-talk modeling) and first-order (with 
cross-talk modeling) domain decomposition method identical. 

5.7. 3D simulations where cross-talk is significant 

In the examples presented at the end of Chapter 4, cross-talk was not included in 
the decomposition methods used (qr-DDM and mk-mask). Through the cross-talk model- 
ing of this Chapter, it is now obvious why cross-talk was not significant in those examples. 
In the first, a 0°/180° alt. PSM with 1:1 dense line/space patterns and 50nm of undercut 
does not exhibit any significant amount of cross-talk between the 0° and 180° phase- wells 
according to Figure 5-9 and Figure 5-10. The same is true for the second example, in which 
a 90°/270° alt. PSM with 1:1.5 semi-dense contacts and 50nm of undercut doesn't have 
high levels of cross-talk as evidenced in Figure 5-10 and Figure 5-12. Consequently the qr- 
DDM and subsequent mk-mask methods were used and yielded accurate image results. 

Two more examples of 2D mask layouts are presented next. The purpose of the first 
example is to demonstrate the successful inclusion of cross-talk terms in 2D layouts, when 
it is erroneous to ignore it, whereas the purpose of the second is to validate once again the 
physical understanding in place for cross-talk, by showing that cross-talk because of the 
corners at the bottoms of the phase-wells is insignificant and can be safely ignored. 

Figure 5-21 shows a 180°/360°, CD=120nm(lX), 1:1 dense layout for contacts 
with no undercut. The rigorous solution of the scattered field at the observation plane for 
the complete alt. PSM under TE (E y ) illumination at A^=193nm is shown in Figure 5-22(a). 
The amplitude of the error incurred using the DDM of Figure 5-2 1(b) is shown in Figure 5- 
22(b). Clearly, the relatively large error of the decomposition method is because cross-talk 
was neglected. Cross-talk effects can be included using only 2D simulations, as, was done 
in Section 5.4. Note that this involves multiple cross-talk terms, both with TE and TM 
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Figure 5-21. Dense contacts: Mask layout and DDM with cross-talk modeling 

(a) 1:1 dense contact hole mask layout with CD=120nm (IX), (b) layout 
decomposition via the DDM into the four constituent single-opening masks 
and (c) schematic diagram of the various cross-talk terms that need to be 
included for an accurate image simulation. 



light, as shown in Figure 5-21 (c): The southwest 180° hole exhibits TE cross-talk effects 
with the southeast 360° hole (and its mirror 360° hole in the x-axis because of periodic 
boundary conditions). The same 180° hole also exhibits TM cross-talk effects with the 
northwest 360° hole (and its mirror 360° hole in the y-axis because of periodic boundary 
conditions). Similar cross-talk effects are observed between the northeast 180° hole and the 
(a) (b) 
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0.96 
pm(4X) 



1.92 



Figure 5-22. Dense contacts: Near field comparison 

(a) Amplitude of the scattered field (Ey) of the mask layout shown in Figure 5-21(a) calcu- 
lated with a rigorous 3D simulation and (b) amplitude of the error of the scattered field of 
the DDM shown in Figure 5-21(b) compared to the rigorous solution in (a). 
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southeast and northwest 360° holes. No cross-coupling occurs between diagonal holes. The 
final aerial images for a system with G=0.3, NA=0.7 and R=4 are shown in Figure 5-23. In 
Figure 5-23(a) the image obtained using the fully rigorous r-mask is shown along with the 
20% intensity contours for the three different approaches (r-mask, DDM without and with 
cross-talk modeling). From the comparison of the three distinct images along the cut-line 
and within the zoom boxes shown in Figure 5-23(b) and (c) it is obvious that including 
cross-talk effects restores the accuracy of the DDM. The normalized mean square error 
within the image field drops by more than a factor of 10, from 23.6E-4 to 1.5e-4. 
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Figure 5-23. Dense contacts: Aerial image comparison 

The image of the mask layout shown in Figure 5-21(a) (X=193nm, a=0.3, NA=0.7, 
R=4) calculated with the rigorous scattered field solution of Figure 5-22(a) is 
shown in (a). Three iso-intensity contours at the 20% level are also shown, one cor- 
responding to the r-mask, one to the qr-DDM without cross-talk modeling and the 
other to the qr-DDM with cross-talk modeling. The three images are compared 
along the indicated cut-line in (b) and within the zoom boxes in (c). 
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In the second example, a 1 80°/360° aJt. PSM layout bearing a pair of square holes 
with CD=120nm (IX) that are diagonally positioned is shown in Figure 5-24(a). The rig- 
orous solution of the scattered field at the observation plane for the complete alt. PSM 
under TE (Ey) illumination at X=193nm is shown in Figure 5-25(a). The amplitude of the 
error if the DDM of Figure 5-24(b) is used is shown in Figure 5-25(b). Observe that the 
level of the error is one order of magnitude smaller than in the previous example. This indi- 
cates that cross-talk effects are insignificant. Although the two holes in the layout are in 
< a > (b) ( C ) 
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Figure 5-24. Diagonal contacts: Mask layout and DDM - insignificant cross-talk 

(a) Mask layout with a diagonally positioned pair of contact holes with CD=120nm 
(IX), (b) layout decomposition via the DDM into the two constituent single-opening 
masks and (c) schematic diagram showing that cross-talk in the diagonal direction 
because of the corners is insignificant. 




0.96 
^m (4X) 




0.96 
urn (4X) 



0.025 



1.92 



Figure 5-25. Diagonal contacts: Near field comparison 

(a) Amplitude of the scattered field (Ey) of the mask layout shown in Figure 5-24(a) calcu- 
lated with a rigorous 3D simulation and (b) amplitude of the error of the scattered field of 
the DDM shown in Figure 5-24(b) compared to the rigorous solution in (a). 
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close proximity, their diagonal displacement puts the sources of cross-talk (edges at bottom 
of phase-wells) in such relative locations, that it is impossible for scattered energy to be 
cross-coupled into neighboring openings. Had the corners been a significant contributor to 
cross-talk effects, the error shown in Figure 5-25(b) would have been comparable to the 
error levels of Figure 5-22(b). The aerial images for a system with a=0.3, NA=0.7 and R=4 
are shown in Figure 5-26. In Figure 5-26(a) the image obtained using the fully rigorous r- 
mask is shown along with the 10% intensity contours for the two approaches (fully rigor- 
ous r-mask and DDM without cross-talk modeling). Comparing the images along the diag- 
onal cut-line and within the zoom boxes shown in Figure 5-26(b) and (c) further reinforces 
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Figure 5-26. Diagonal contacts: Aerial image comparison 

The image of the mask layout shown in Figure 5-24(a) (A=193nm, a=0.3, NA=0.7, 
R=4) calculated with the rigorous scattered field solution of Figure 5-25(a) is shown in 
(a). Two iso-intensity contours at the 10% level are also shown, one corresponding to 
the r-mask and the other to the qr-DDM without cross-talk modeling. The two images 
are compared along the indicated cut-line in (b) and within the zoom boxes in (c). 
Cross-talk effects are negligible, regardless of the fact that the contacts are in close 
proximity along the diagonal direction. 
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the idea that corner cross-talk effects are negligible. The two images are hardly distinguish- 
able since the normalized mean square error of the DDM without cross-talk modeling is 
6.2E-7! 

5.8. Off-axis illumination 

Here, it is of great importance to check to validity of the cross-talk models for off- 
axis illumination that would be required for the application of the decomposition methods 
in imaging formulations such as Abbe's or a modified Hopkins' method with varying dif- 
fraction orders vs. angle of incidence. 

In the examples of Figure 5-27 a 180°/360°, 1:1 dense line/space alt. PSM with 
CD=400nm (4X) and vertical or tilted sidewalls (same geometries as in Figure 5-6) is illu- 
minated by a plane wave traveling at -10.9° off-normal. Omitting the cross-talk effects 
(O^-order DDM) is inaccurate, as shown in Figure 5-27 (a) arid (b). Including the cross-talk 
terms shown in Figure 5-27(c) and (d) leads to an accurate result for the mask diffraction 

(l st -order DDM) that is almost indistinguishable from the rigorous r-mask model (complete 
alt. PSM geometry). 
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Figure 5-27. Examples of l st -order DDM with off-axis illumination 

The examples of Figure 5-6 are revisited to demonstrate the validity of the cross-talk 
model in off-axis illumination. The illuminating plane wave strikes the mask at 10.9° 
off the normal (z-axis). A comparison of the three simulation approaches is shown: 
Complete alt. PSM simulation (no decomposition, r-mask), O^-order and l st -order 
DDM for vertical (a) and sloped (b) sidewall profile. The normalized amplitude of CT 12 
and CT 2 i is also shown for vertical (c) and sloped (d) sidewall profiles. 

5.9. Conclusions 

In this Chapter a novel methodology to model cross-talk between phase-wells in 
alt. PSMs was presented. Through the calculation of the scattered field from an isolated 90° 
air/glass discontinuity a new model for cross-talk was proposed. According to this model, 
the edges (not the faces or the corners) at the bottom of the phase-wells are the primary 
sources of cross-talk. The validity of the model was established directly through the suc- 
cessful first-order decomposition of alt. PSM in both 2D and 3D simulations that yield 
accurate results when compared with the simulation of the exact geometry. Also, indirect 
evidence of the validity was seen in the simulation of the half-covered sidewalls of the 
phase-wells and the M0O3 ARC of Section 5.6. Cross-talk can potentially be a source of 
image imbalance but can also be used to "correct" imbalance. Recipes that reduce cross- 
talk and the way they achieve this task were also presented. Both these ideas are manufac- 
turable (specifically the Cr-covered sidewalls has been demonstrated by Levenson et al. 
[50]). It was cautioned however that the elimination of cross-talk although gets rid of a 
hard-to-control phenomenon does not necessarily lead to superior images and cross-talk is 
definitely not the sole contributor to image imbalance. Although light power owing to 

121 



c c 

cross-talk is only a small portion of the total signal power (-1% or less), the abrupt phase 
changes in a PSM can amplify the effects of cross-talk through the constructive or destruc- 
tive interference of cross-talk with the near diffracted field. The cross-talk model can also 
handle the important case of off-axis illumination. 
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The Edge Domain 




A powerful extension of the domain decomposition methods developed in Chapter 
4 for the rapid and accurate simulation of light scattering from advanced photomasks is 
demonstrated here. 

Researchers have in the past investigated the decomposition of mask layouts into 
individual openings, and even into edges, for ideal thin masks under coherent illumination 
[65], but such techniques have never been applied to the rigorous solution for the scattered 
field. Domain decomposition techniques, where a large electromagnetic problem is broken 
up in smaller pieces and the final solution is arrived at by synthesizing (field-stitching) the 
elemental solutions, have been recently proposed for the study of one-dimensional binary 
(phase only) diffractive optical elements [49], [72]. Another research group working on the 
same problem demonstrated how to create and use a perturbation model for binary edge- 
transitions based on the product of the ideal, sharp transition and the continuous field vari- 
ations in the vicinity of the edge [48]. Independently, a similar technique to field-stitching 
has been invented and presented in Chapter 4 and Chapter 5 for the simulation of two- 
dimensional layouts of advanced photomasks (alternating PSM, masks with OPC) that 
properly models interactions from neighboring apertures (cross-talk) and furthermore 
takes advantage of the spectral properties of the diffracted fields to come up with a compact 
model for the edge-transitions [2]. 

In this Chapter, a domain decomposition method based on edges is described, 
where the pre-calculated electromagnetic field from the diffraction of isolated edges is 
recycled in the synthesis of the near diffracted field from arbitrary two-dimensional dif- 
fracting geometries. Although at first glance it might be tempting to discard edge-decom- 
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position techniques as inaccurate for masks with large vertical topography and features 
with small lateral dimensions, it will be shown that for a large set of practical situations they 
are accurate. A key difference in this method from the one in [48] is that the sum of the 
complex fields is taken instead of the product. The reason for that lies in the linearity of the 
Kirchhoff-Fresnel diffraction integral [15] and it is exactly because of this linearity that the 
decomposition method is possible. 

The Chapter is organized as follows: First, a new method for the mask domain 
decomposition is introduced and applied in simple ID layouts. The limits of this method 
are systematically examined and presented next. Then, the method is expanded in 2D lay- 
outs and its algorithmic implementation for arbitrary layouts is considered. Simple 2D lay- 
outs that establish the validity and accuracy of the approach and large, arbitrary layouts that 
demonstrate its versatility are covered next. Finally, an important direction for further 
enhancing the efficiency of the method that is based on spectral matching is discussed. 

6.1. The edge-DDM applied in ID layouts 

A natural extension of the DDM shown in Figure 4-3(b) that is also based on the 
linearity of the Kirchhoff-Fresnel diffraction integral is shown in Figure 6-1. The ID mask 
layout is first decomposed into single opening masks as in Figure 4-3(b) and subsequently 
each opening is decomposed into two edges and a uniform illumination field as in Figure 6- 
1 . Such a decomposition can be proven to be exact for infinitely thin and perfectly conduct- 
ing screens through the use of the electromagnetic form of Babinet's principle [16], but it 
should break down if applied to photomasks with "thick" vertical structures (compared to 
A.) and not perfectly absorbing screen materials. Also, because of the reduced ability of light 
to penetrate through small (compared to X) openings it is expected to break down for small 
mask features. This decomposition method will be hereafter referred to as the edge domain 
decomposition method or edge-DDM. The application of the edge-DDM to a single isolated 
space and an isolated line are shown next, followed by a systematic evaluation of the limits 
of applicability of the method. 
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Figure 6-1. Decomposition of an opening into two edges 

An opening with transmission T and phase <(> in an infinitely thin, opaque screen 
can be decomposed by way of the linearity of the Kirchhoff-Fresnel diffraction 
integral into two edges. The uniform plane wave illumination is subtracted to 
restore the light level everywhere. 

6.1.1. Decomposition of a space into two edges 

Consider the isolated mask space (opening) shown in Figure 6-2(a). The amplitude 
of the complex E y -field under E y (TE) illumination is shown in (c) and the scattered field 
across the observation plane is shown in (e). According to Figure 6-1, the field across the 
observation plane can be alternatively obtained by using the scattered field from an edge 
that has the same vertical profile as the space, as shown in Figure 6-2(b) and (d). The result- 
ing field from this edge-DDM is overlaid on the plot of Figure 6-2(e) along with the field 

from the edge-scattering. The normalized mean square error 1 (NMSE) between the r-mask 
and the edge-DDM is 0.19%. In Figure 6-2(f) the spectra of the r-mask and edge-DDM are 
compared and it is seen that they differ by 0.13% in a NMSE sense if the whole spectrum 
of propagating waves is considered and by merely 0.007% within the collection ability of 



1. The normalized mean square error is defined by: NMSE = J|£ r - E^dS/\\E^[dS , where E r is 

the "true" field resulting from a rigorous simulation of the complete mask geometry and E s is the syn- 
thesized field via the edge-DDM. 
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a projection system with NA=0.8 and R=4. Since only the portion of the spectrum that is 
collected by the projection system contributes to the image formation, the exceptional accu- 
racy of the edge-DDM is this regime is more than adequate for accurate image simulations. 

6.1.2. Decomposition of a line into two edges 

Similar steps can be followed for the isolated line shown in Figure 6-3(a) illumi- 
nated with E x (TM) polarized light. Note that it is not typical to have a line as shown here, 
where in both sides the glass has been etched. Normally only one side of the line is phase- 
shifted, but this example simply aims to introduce the principle of the edge-DDM. A line 



( a > (b) 




^1 (4X) 1^ ( nm) -i 

Figure 6-2. Edge-DDM applied on an isolated space 

Geometry of mask with isolated space (a) and edge (b). Note that the vertical 
profile of the edge is the same as that of the isolated space. Amplitude of com- 
plex field from rigorous TEMPEST simulation for space (c) and edge (d). The 
true scattered field (r-mask) across the observation plane in (c) is compared with 
the edge-DDM (synthesized field) in (e) and their spectra are compared in (f). 
The normalized mean square errors (NMSE) are indicated on the plots. 
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with one side phase-shifted and the other left intact can still be decomposed via the edge- 
DDM using the diffraction fields from the two different edges. The amplitude of the com- 
plex E x field everywhere in the domain of Figure 6-3(a) is shown in (c) and for the edge of 
Figure 6-3(b) in (d) respectively. The true (r-mask) scattered field is compared with the 
synthesized field with edge-DDM in (e) and the spectra are compared in (f). The error 
(NMSE) is seen to be small in both the near field and spectra plots and specifically the error 
of the through-the-lens (TTL) spectrum is just 1.5E-5. 
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Figure 6-3. Edge-DDM applied on an isolated line 

Geometry of mask with isolated line (a) and edge (b). Note that the vertical pro- 
file of the edge is the same as that of the isolated line. Amplitude of complex 
field from rigorous TEMPEST simulation for line (c) and edge (d). The true scat- 
tered field (r-mask) across the observation plane in (c) is compared with the 
edge-DDM (synthesized field) in (e) and their spectra are compared in (f). The 
normalized mean square errors (NMSE) are indicated on the plots. 
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6.1.3. Limits of the edge-DDM 

In the examples of Figure 6-2 and Figure 6-3 the edge-DDM is seen to be accurate 
compared with the simulation of the exact mask structure. Although the vertical mask 
topography (180° - 170/im = IX) is comparable to the wavelength, the fact that the lateral 
dimensions of the space and the line are relatively large compared with the wavelength 
(CD = 400/1/71 > 2 x 193/im = 2 x X ) renders the edge-DDM accurate. One expects that 
as the lateral sizes become smaller and the vertical mask topography even larger the edge- 
DDM will eventually break down. The limits of the edge-DDM were investigated through 
exhaustive simulation of isolated spaces and lines on the mask with various etched depths 
and amounts of undercut. The sidewall aneles were alwavs k-pnt 

at the bottom of the etched wells was included. The mask CD was varied in the set {2jtm, 
lum, 0.6um, 0.4um, 0.2fim(~U), O.lum}, the etched depth in {0°, 90°, 180°, 270°, 360°} 
and the amount of underetch in {Onm, 25nm, 50nm, 75nm, lOOnm}. All 2 (space/line) by 
2 (TE/TM illumination) by 6 (CD sizes) by 5 (etched depths) by 5 (underetch) = 600 cases 
were run and the rigorous (r-mask) scattered field was automatically compared with that 
obtained from edge-DDM. The results are spectacular! The edge-DDM can readily achieve 
better than 1% accuracy (in a NMSE sense) in the near field for CD=400nm (4X) or larger 
regardless of depth, underetch and polarization. The accuracy in the TTL-spectrum is suf- 
ficient in those cases even for inspection (R=l) simulations. For lines, the edge-DDM 
maintains excellent accuracy in the TTL-spectrum with R=4 down to a mask CD of lOOnm 
(-0.5A) with 90° of etched depth and Onm of undercut, or mask CD of 200nm (~U) with 
270° of etched depth and 25nm of undercut. For spaces, the accuracy of the method is 

acceptable at least down to mask CD of 200nm (~1A) with 270° of etched depth and 50nm 
of underetch. 

Current state-of-the-art photomasks and processes utilize sub-resolution assist fea- 
tures that are never less than one wavelength in size on the mask. This clearly implies that 
the edge-DDM has the potential to accurately simulate the most advanced photomask tech- 
nologies (OPC, alt. PSM, OPC on alt. PSM). Again, as with the DDM of Chapter 4, the true 
power of the method lies on its application in 2D mask layouts. However, one key differ- 
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ence with the DDM of Chapter 4 is the following: Only one rigorous 2D simulation of an 
edge suffices for the scattered field reconstruction of any arbitrary size on the mask! The 
DDM of Chapter 4 requires a separate simulation of every different mask size. Next, the 
application of the edge-DDM in 2D layouts is developed. 

6.2. The edge-DDM applied in 2D layouts 

The application of the edge-DDM in 2D layouts is straightforward. A mask layout 
comprises of a large number of edges positioned at different locations and having different 
orientations. No matter how complicated the mask technology, there are usually only a 
small number of different types of edges present in the layout. For example, a single expo- 
sure 0°/90°/270° alt. PSM has five types of edges: i) Cr-layer/0° edge, ii) Cr-layer/90° 

edge, iii) Cr-layer/270° edge, iv) 0°/90° edge and 0°/270° edge. Depending on the orienta- 
tion of each edge in the layout it "sees" and responds to the incident field differently. If the 
illuminating field is a TE (Ey) normally incident plane wave and the edge is oriented along 
the y-axis it "sees" TE illumination, but if it is oriented along the x-axis it "sees" TM illu- 
mination. This is illustrated in Figure 6-4. The response of each edge to its respective illu- 
mination is taken into account rigorously from a pre-stored 2D edge-diffraction simulation. 
However, corner effects are effectively ignored, since the finite extent of each edge is not 
rigorously taken into account. Instead, the scattered field at the end points of every edge is 
abruptly terminated (truncated) in a perfect square-wave fashion to the field value of the k- 
mask model. As seen in Section 4.3.4, corner effects in typical imaging situations are 
mapped at the extremities of the spectra and do not contribute to the image formation. 
Owing to this observation, the edge-DDM applied in 2D layouts is successful in accurately 
capturing the true electromagnetic behavior of edges in a rapid manner, while the insignif- 
icant corner effects are safely ignored. 
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Figure 6-4. Edge-DDM applied on an arbitrary Manhattan-type polygon 



Depending on the orientation of each edge in the layout it "sees" and responds to the 
incident field differently. Here the illuminating field is a TE (E y ) normally incident 
piane wave and the edges that are oriented along the y-axis "see" TE illumination 
(parallel to the edge), whereas the edges that are oriented along the x-axis "see" TM 
illumination (perpendicular to the edge). A separate 2D edge simulation is required to 
capture the different response of the edge to different field polarizations. 

Slightly different embodiments of the edge-DDM applied in 2D layouts can be 
devised. One that is powerful and general in its application is described next. Then, scatter- 
ing results from various edge profiles that would be encountered in alt. PSMs are presented 
for completeness. Next, application of the edge-DDM in simple 2D mask layouts is pre- 
sented and the results of the method are compared with fully rigorous r-mask models in 
order to establish the validity of this approach. Finally, a large portion of a real layout from 
a single exposure 0°/90°/270° alt. PSM that is too big for 3D rigorous simulation is simu- 
lated via the edge-DDM and the aerial image is compared with the simple k-mask model. 

6.2.1. Algorithmic implementation of the edge-DDM 

First, every distinct edge that is present on the layout is pre-simulated with all pos- 
sible illumination directions (field polarizations) that are required, based on the orientations 
of the edge encountered within the layout. For simplicity only Manhattan-type layouts are 
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considered, where the edges are oriented along either the x- or the y-axis and require elec- 
tromagnetic simulation under both parallel and perpendicular to the edge polarizations. 
The extension to non-Manhattan layouts although much more computationally laborious 
is straightforward. Subsequently, each distinct mask layer is broken up into a set of mutu- 
ally disjoint (non-overlapping) rectangles. For example, a plain binary mask consists of 
two layers, namely the absorption layer (typically Cr-based) and the clear (uncovered) 
layer, whereas a more advanced 0°/90 o /270° three-phase alt. PSM consists of four layers, 

namely the absorption layer and the 0° clear, 90° clear and 270° clear layers. Next, looping 
through all rectangles of each layer, the type of each one of the four edges of the rectangle 
is determined based on what its neighboring layer is. Trivial edges (when the neighboring 
rectangle is of the same type) are discarded. It is also possible that only part of a single edge 
from a rectangle is neighboring with one layer and the rest is neighboring with a different 
layer (or different layers). In that case the edge is broken up in multiple edges, in a way that 
tracks the neighboring layer. What is left at the end of this process is a set of all non-trivial 
edges that are present in the layout (location, size and orientation) and their type, i.e. which 
are the two mask layers on each side of the edge. Given the mask layout, the k-mask 
(Kirchhoff-mask) model is readily available. Finally, the respective (complex) difference 
of the true edge-scattering from the ideal and sharp k-mask model is added to all non-trivial 
edges accordingly. This revised mask model that results from the edge-DDM includes 
accurate information about the electromagnetic scattering from the edges and it is a quasi- 
rigorous mask model (qr-mask) in the same sense that the term was used in Section 4.3, 
where only 2D mask simulations were used to approximate a 3D electromagnetic problem. 

The above algorithm was implemented in the MATLAB environment. The part that 
deals with the polygon and edge extraction can probably be judged as rudimentary and 
most likely cannot compete in speed with more sophisticated implementations in produc- 
tion caliber CAD software, but it suffices for the purposes of this work, that is, to demon- 
strate the proof of concept of the edge-DDM and the feasibility of including accurate edge- 
scattering information in rapid aerial image simulations. 
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6.2.2. Dependence of edge-scattering on profile and polarization 

Examining the departure of the scattered field from an edge from its ideal, sharp 
edge-transition and how it depends on the profile of the edge and the incident field polar- 
ization can provide intuition of the physical mechanisms involved. 

Figure 6-5 shows the amplitude of the (near) scattered field from an edge with var- 
ious profiles. The ideal, sharp edge-transition is also shown on every graph for comparison. 
In Figure 6-5(a) and (b) the edge-scattering for different etched depths (0°, 180° and 360°) 
and vertical sidewalls with no undercut is shown for TE (Ey) and TM (E x ) illumination 
respectively. The effect of undercut is depicted in Figure 6-5(c) and (d), where the edge- 
scattering for different amounts of undercut (Onm, 25nm, 50nrn and 75nm) and vertical 
sidewalls with 180° etched depth is shown for TE and TM illumination respectively. It is 
interesting to note that increasing the etched depth is "making the edge response slower", 
meaning that the near scattered field from the edge requires more distance to reach its "on" 
value from its "off* value, although this effect is more pronounced for the first increment 
from 0° to 180° rather than from 180° to 360°. Similarly, increasing the amount of undercut 
appears to "make the edge response faster" and furthermore this behavior appears to be lin- 
ear, since for every increment of the undercut the edge-transition gains a constant distance 
in its "off-on" behavior. 

It is instructive to isolate the difference of the true edge-scattering from the ideal, 
sharp edge-transition. This is done in Figure 6-6, which depicts the amplitude of this (com- 
plex) difference for all respective plots of Figure 6-5. Observe that, indeed, by increasing 
the etched depth the difference has a larger spatial extent, whereas increasing the amount 
of undercut confines the spatial extent of the difference. This is crucial, because it is exactly 
the spatial extent of the difference that affects its spectral distribution. If the difference is 
well confined, its spectral distribution is spread out to higher spatial frequencies and the 
lower frequency content is small. In such a case the complex difference can be neglected 
and the approximation of the ideal, sharp transition is adequate. 
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Figure 6-5. Dependence of edge-scattering on profile and polarization 

Normalized amplitude of the near field across the observation plane for an edge with 
0°, 180° and 360° of etched depth under TE (a) and TM (b) polarizations. Normalized 
amplitude of the near field for an edge with 180° of etched depth and Onm, 25nm 
50nm and 75nm of undercut under TE (c) and TM (d) polarizations. The ideal sharp 
edge-transition is shown in all cases for comparison. Observe that the "edge 
response" becomes slower for increasing etched depths and faster for increasing 
undercuts. 
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Figure 6-6. Complex difference of true edge-scattering from ideal, sharp transition 

The amplitude of the complex difference of the true electromagnetic field scattered 
from an edge from the ideal, sharp transition (step-function) is shown in (a)-(d) for all 
edge profiles and illumination polarizations of Figure 6-5 respectively, (a) and (c) are 
for TE polarization and (b) and (d) are for TM polarization. The spatial extent of the 
difference grows with increasing etched depths, as seen in (a) and (b), and shrinks 
with increasing undercuts, as seen in (d) for TM polarization. The behavior in (c) 
appears somewhat erratic, but if the phase is also taken into account, it is shown in 
Figure 6-7 that this is also the case for TE polarization. 



However, when the difference has a larger spatial extent, then its lower frequency 
content becomes significant and it cannot be neglected. This is shown in Figure 6-7, which 
depicts the magnitude of the spectra of all respective difference distributions of Figure 6-6. 
Observe that the TTL spectrum of the difference attains larger values when the etched depth 
increases and smaller values when the amount of undercut increases, as expected from 
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Figure 6-5 and Figure 6-6. Remember that the normalization used in spectra plots is that a 
constant (DC) value of 1 (clear field) would produce a Dirac 5-function at zero frequency 

with amplitude of 19.4 (Jr\ = Jj\i 0 /z Q ~ 19.4Q). 




(c) (d) 




Figure 6-7. Spectrum of complex difference 



The spectra of the complex differences shown in Figure 6-6 are depicted here respec- 
tively, (a) and (c) are for TE polarization and (b) and (d) are for TM polarization. 
Observe that the through-the-lens (TIL) parts of these spectra, which for typical 
lithographic imaging situations (R=4,NA~0.7-0.8) extent up to - +/- 0.008 (nm)" 1 are 
smaller for smaller etched depths and larger undercuts, regardless of polarization, as 
expected from the observations on the plots of Figure 6-5 and Figure 6-6. 
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6.23. Application of edge-DDM in simple 2D layouts 

The application of the edge-DDM on an isolated square hole and an isolated square 
island (post) and comparison of the results with fully rigorous 3D mask simulations (r- 
masks) are presented here. The edge profile has a glass etched depth of 180° and an under- 
cut of 50nm. The dimensions of both the hole and the island are 400nm x 400nm (4X), cor- 
responding to approximately 2x2 wavelengths, for X=193nm. 

Figure 6-8(a) and (b) depict the amplitude of the scattered field across the observa- 
tion plane below the hole obtained via the edge-DDM and via fully rigorous 3D simulation, 
respectively. The amplitude of their difference (error of the edge-DDM) is shown in 
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Figure 6-8. Edge-DDM applied on an isolated square hole (scattered field) 

(a) Synthesized with the edge-DDM and (b) rigorously calculated scattered field 
(amplitude), across the observation plane below a 400nm x 400nm isolated hole. A 
180° etch and 50nm undercut has been applied. The amplitude of the error (difference 
of r-mask from edge-DDM) is shown in (c). Locally the amplitude of the error 
reaches the 50% level, but the normalized mean square error does not exceed 5%. 
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almost 50% of the clear field value. The normalized mean square error is however more 
contained and is only about 4.9% if calculated within a 3|um by 3^im area or 3.6% if cal- 
culated within the 0.4p,m by 0.4jim open area of the hole, as shown in Figure 6-8. From the 
work presented so far, it should be obvious that neither of these error levels is immediately 
relevant in the simulation of the image formation. What is rather more relevant is the error 
incurred on the diffraction orders that are collected by the imaging optics. If the spatial fre- 
quencies of the error depicted in Figure 6-8(c) are sufficiently high, this error will be dis- 
carded by the projection lens. This concept is shown in Figure 6-9, where the spectra 
(magnitude) of the scattered fields from edge-DDM and from rigorous 3D simulation are 
shown in (a) and (b) respectively. The marginal circles indicating the propagating plane 

waves (with \k\ <2n/X) and the collected by the imaging system plane waves (with 



(a) Propagated orders (b) 




Mnm)-» 

Figure 6-9. Edge-DDM applied on an isolated square hole (spectrum) 

Magnitude of the plane wave spectra of the near field plots in Figure 6-8. Spectrum of 
the synthesized (edge-DDM) field in (a) and of the rigorous field (from 3D simula- 
tion) in (b). The spectrum of the error is depicted in (c). The error level within the col- 
lected orders circle is seen to be an insignificant 0.02%. 
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|l| < 2nNA(l + a)/XR), for NA=0.75, R=4 and a=0.3, are indicated on these plots. The 
magnitude of the error in the spectrum of the field obtained from edge-DDM is shown in 
Figure 6-9(c). Observe that most of the error is concentrated in the higher spatial frequen- 
cies for |k x | approaching 0.03 and only small levels of error exist in the collected orders. 
The normalized mean square errors of all propagated and collected orders are 2.2% and 2E- 
4 respectively. 

Similarly, the amplitude of the scattered fields across the observation plane below 
the isolated island obtained with the edge-DDM and with fully rigorous 3D simulation and 
the amplitude of the error are shown in Figure 6-10. Again, a locally significant error 
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Figure 6-10. Edge-DDM applied on an isolated square island (scattered field) 

(a) Synthesized with the edge-DDM and (b) rigorously calculated scattered field 
(amplitude), across the observation plane below a 400nm x 400nm isolated island 
(post). A 180° etch and 50nm undercut has been applied. The amplitude of the error 
(difference of r-mask from edge-DDM) is shown in (c). Locally the amplitude of the 
error reaches the 50% level, but the normalized mean square error does not exceed 
2%. 
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amplitude near the corners is evident. However, the energy of this error is concentrated 
mostly at the higher spatial frequencies, for |k x | approaching 0.03 and the error in the col- 
lected orders is seen to be small, only 0. 16% in a NMSE sense. Note that the large DC com- 
ponent of the isolated island mask layout was subtracted in order to make meaningful 
comparisons of the edge-DDM with fully rigorous simulations. Otherwise this large DC 
component would artificially show an even better accuracy, since the error level in com- 
parison with a much stronger signal (containing the DC term) would be smaller. 
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Figure 6-11. Edge-DDM applied on an isolated square island (spectrum) 

Magnitude of the plane wave spectra of the near field plots in Figure 6-10. Spectrum 
of the synthesized (edge-DDM) field in (a) and of the rigorous field (from 3D simula- 
tion) in (b). The spectrum of the error is depicted in (c). The error level within the col- 
lected orders circle is only 0.16%. 



It is worth mentioning that each fully rigorous mask simulation of the patterns of 
Figure 6-8 and Figure 6-10 that was performed with TEMPEST used approximately 
1.8Gbytes of memory and required almost 2 days on a 550Mhz CPU. On the other hand, 
once all of the required 2D edge simulation results are available the scattered field can be 
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obtained with the edge-DDM in less than 1 sec! The simulation of the isolated edge-diffrac- 
tions is an off-line procedure. It only needs to be performed once and the edge-diffraction 
results are recycled as many times as are necessary for the reconstruction of the scattered 
fields of arbitrary layouts. In the above examples, it took approximately 5min for each rig- 
orous, isolated-edge 2D simulation. 

6.2.4. Example of edge-DDM on a large, arbitrary layout of a 0°/90°/270° alt. PSM 

The edge-DDM is used here to generate the scattered field from a three level (0°/ 
90°/270°) alt. PSM. A portion of the arbitrary target layout is shown in Figure 6-12 before 
and after the phase assignment, The rectangle and edge extraction according to the algo- 
rithm presented in Section 6.2.1 is performed within the 3^im by 4fim (IX) coherence 
window shown and the scattered field obtained with the edge-DDM inside the 12}xm by 




Figure 6-12. Example: Edge-DDM on a large 2D mask layout (algorithm and near fields) 

The target mask layer and the mask after assignment of an arbitrary layout are shown. A 
three phase-level alt. PSM is assumed (0°/90°/270°). The rectangle and edge extraction 
within the 3\im by 4|Lim (IX) coherence window for the various polygons is also shown. 
Finally, the scattered near field within the coherence window for both TE and TM 
(X=I93nm) illumination is depicted. 
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16pim (4X) coherence window is also shown in Figure 6-12 for TE and TM illuminations 
at A=193nm. Note that a fully rigorous 3D mask simulation of a 12^m by I6\im (4X) 
layout would require close to 40Gbytes of memory and would only be feasible on a multi- 
CPU architecture. The scattered field is obtained with the edge-DDM implementation in 
less than lmin. This time includes the polygon and edge extraction and the field synthesis 
from pre-calculated edge-diffractions. Therefore, the results of the edge-DDM cannot be 
immediately compared with a fully rigorous solution, but the accuracy of the edge-DDM 
within the imaging system's collection ability can reveal the inaccuracies of the ideal k- 
mask, where all transitions are assumed perfectly sharp. This is done in Figure 6-13, where 
the images (at the 30% intensity level) obtained with the k-mask and the edge-DDM for 
two different undercuts (Onm and 50nm) are compared. The system parameters are 
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Figure 6-13. Example: Edge-DDM on a large 2D mask layout (images) 

Image (30% intensity contours) of the layout within the 3jim by 4jiim (IX) coherence 
window of Figure 6-12. The ideal k-mask (Kirchhoff-mask) model is compared with 
the edge-DDM for an undercut of Onm and 50nm. The edge-DDM is capable of 
revealing the inaccuracies of the k-mask model and of accurately taking into account 
the edge-scattering of different edge profiles. 
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A=193nm, NA=0.7, a=0.3 and R=4. The inaccuracies of the k-mask mode] are revealed at 
various locations in the image as seen in Figure 6-13. The effect on the image of changing 
the undercut from Onm to 50nm can be quickly and accurately evaluated. Observe the 
subtle differences evident in the intensity contours of the Onm and 50nm masks of Figure 6- 
13. 

6.3. Further speeding-up the edge-DDM 

Although the edge-DDM presented earlier achieves excellent accuracy and tremen- 
dous speed-ups as compared to a fully rigorous mask simulation, there is still one more 
parameter to exploit in order to further speed-up the method. As was shown in Chapter 4, 
uiv iiiai^nwu x^ii^iniuii maajv niuuci ^niK-iiiiisK; was ueveiopea, mere is no need to 
carry all the details of the true field scattering as far as the imaging process is concerned 
(with larger than 1 reduction factors), since the higher spatial frequencies are not used. If a 
lumped parameter model can be devised that closely matches the TTL spectrum of the scat- 
tered field then the task of an accurate field representation is accomplished. In Chapter 4 a 
rect-function, a raised cosine function and a Gaussian function were seen to successfully 
mimic the TTL spectra of scattered fields from various size and profile openings. A similar 
principle can be applied here with the purpose of matching the TTL spectra of edge-dif- 
fracted fields. Assuming a restriction to piecewise constant models, the scattered field from 
an isolated edge under TE or TM illumination can be approximated using the multi-step 
functions shown in Figure 6-14. Although not shown in Figure 6-14, the amplitude, phase 
and size of all steps are optimized to achieve the best spectral matching of the collected 
(TTL) orders with the orders resulting from the continuous edge-diffracted field. Note that 
the key factor for a good spectral match is to capture the finite extent of the "off-on" edge- 
transition. A piecewise constant function is probably not the best choice to achieve this 
task, but on the other hand a piecewise constant function offers a good lumped parameter 
approach and it lends itself nicely as an extension of the unmatched Kirchhoff approach, 
where all edge-transitions are assumed perfectly sharp. 

If the matched bandwidth edge-DDM is applied on the 400nm by 400nm (4X) iso- 
lated hole of Figure 6-8 the scattered field across the observation plane can be represented 
in the compact way shown in Figure 6-15(a). Recall that the edges parallel to the y-axis 
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Figure 6-14. Matched bandwidth approximation to edge-diffraction 

Instead of approximating the true edge-diffraction by a single step function (geomet- 
rical or Kirchhoff approximation) a multi-step piecewise constant function is 
employed. The amplitudes, phases and sizes of all steps are optimized to achieve the 
best spectral matching of the collected orders with the orders resulting from the con- 
tinuous edge-diffracted field. The matched bandwidth approximations to the field dif- 
fracted by a 180° edge with 50nm undercut under TE and TM illuminations are 
shown in (a) and (b) respectively. j 

"see" TE polarization and the edges parallel to the x-axis TM polarization. The multi-step 
matched BW approximations to the edge-diffraction shown in Figure 6-14 for TE and TM 
illuminations intermingle in the x- and y-directions to produce the multi-color mask of 
Figure 6-1 5(a), where each color represents a different amplitude and phase region of the 
mask. The spectrum of the error of this approach is depicted in Figure 6- 15(b). Note that 
compared to the error shown in Figure 6-9(c) this error has increased by almost four-fold. 
Nevertheless, the error in the TTL spectrum still remains at low enough levels for accurate 
imaging simulations. A better approximating function, other than the piecewise constant 
step function, that is more capable of capturing the true edge-diffraction profile (such as 
linear or higher order polynomial) should result in better spectral matching and conse- 
quently more accurate imaging simulations. 



143 



c 



(a) 



(b) 




NMSE(propagated orders) = 5.7% (was 2.2%) 
NMS£(collected orders) = 8.6e-4 (was 2e-4) 

Figure 6-15. The matched BW edge-DDM applied to the hole of Figure 6-8 

(a) Representation of the scattered field below a 400nm by 400nm hole (etched depth 
180°, undercut=50nm) using the matched bandwidth edge-DDM, as shown in 
Figure 6-14. Each color represents a different amplitude and phase region of the 
mask. The spectrum of the error of the representation shown in (a) compared to the 
true solution of Figure 6-8(b) is depicted in (b). The error level within the collected 
orders circle is seen to increase from 0.02% to 0.086%. 

6.4. Conclusions 

This Chapter has dealt with the development of the edge-DDM (edge domain 
decomposition method). This method is an extension of the domain decomposition meth- 
ods proposed and developed in Chapter 4 that introduces tremendous versatility, since only 
a small number of isolated edge-diffraction simulations is shown to contain all the neces- 
sary information for the synthesis of the scattered field from arbitrary 2D mask layouts and 
subsequent accurate imaging simulations. The limits of this method are reached when the 
mask features are smaller than a wavelength in size and the vertical mask topography is 
large. Through a systematic process it was determined that features as small as at least 
200nm bearing 270° of phase-wells can be accurately decomposed via the edge-DDM, at 
A,=193nm. The method was tested for simple 2D layouts, where rigorous mask simulations 
are possible. Excellent accuracy accompanied by speed-up factors of 172,800 (lsec vs. 
2days) were demonstrated. The accuracy of the method was attributed to the fact that on 
one hand the edge-diffraction phenomena are modeled rigorously through 2D edge-diffrac- 
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tion simulations and that on the other hand the errors incurred-during the synthesis (prima- 
rily near the corners) are mapped at the extremities of the spectrum of propagating plane 
waves and do not contribute to the image formation. It was emphasized that all necessary 
rigorous 2D simulations of the diffraction from isolated edge profiles is performed off-line 
and recycled for the diffraction calculation of arbitrary layouts. If needed, the accuracy of 
the edge profile diffraction simulations can be pushed to extreme limits, since it is per- 
formed only once. This can aid for example in the correct simulation of unusual edge pro- 
files whose geometrical details would require excessive discretization of the domain. It 
was also shown that the exact details of the edge-diffraction are not necessary for accurate 
imaging simulations (with larger than 1 reduction factors) and a piecewise constant, multi- 
step edge-transition model that matches the TTL spectrum was introduced. The range of 
validity of the edge-DDM is expected to be appropriate for rapid and accurate evaluation 
of aerial images whenever speed is critical, as for example in full-chip OPC software and 
die-to-database comparisons in the inspection of alt. PSMs. 
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7 Characterization of Phase Defects in 
Phase Shift Masks 



PHASE-SHIFTING MASKS (PSMs) have shown promise in increasing the resolution 
and depth of focus in photolithography. Alternating phase-shifting masks (alt. PSMs) are 
one of the two most widely accepted types of PSM (the other being attenuated). Their fab- 
rication process involves at least one quartz etching step, during which the patterned binary 
(Chromium on quartz) mask attains the required 180° phase difference between alternating 
sides of the Chromium-covered quartz. During the quartz etching process phase defects can 
occur, which are generally arbitrarily shaped three dimensional posts of quartz or voids in 
the clear area. 

The strange electromagnetic behavior of light passing through phase-edges exam- 
ined in Chapter 5, where it tends to emanate from the high refractive index material and 
deplete the low refractive index material ([105], [109]) can aid in understanding the impact 
of phase defects, since at a first glance they can be viewed as a number of phase-edges in 
close proximity. Recent simulation and experimental studies have related the deterioration 
of the feature linewidth with the size of phase defects and have also given an example of 
the fundamental difference between the two types of phase defects, namely quartz posts and 
voids in the quartz [23], [93], [30], [33], [4]. 

The main two themes of this Chapter are the following: First, imaging simulations 
that utilize the fully rigorous r-mask model are used to assess the tendency of non-planar 
phase defects to print in optical lithography with alternating phase-shifting masks. A sche- 
matic diagram of a defective alt. PSM is shown in Figure 7-1, where both types of defects 
are depicted. The steady state electromagnetic field throughout the volume of the mask 
structure is calculated with TEMPEST under the excitation of a harmonic, linearly polar- 

146 



c 



c 



Uluimnation plane_ 



glass voi d 



quartz 

180° phase edge 



observation plane (near field) 
vacuum (air) 




Figure 7-1. Schematic diagram of a defective alt. PSM 

Both types of phase defects are shown (voids and posts). The 180° phase edge is 
located at the center behind the Chromium feature (50% underetch). Incident light 
originates from the illumination plane and the full vector electromagnetic solution is 
calculated with TEMPEST everywhere in the domain. 

ized and normally incident plane wave, with illumination wavelength ^=248nm, originat- 
ing from the illumination plane. The diffraction orders at the observation plane are then 
extracted as described previously and form the input to the imaging program (SPLAT) in 
order to calculate the aerial image. The imaging system assumed in this first part is aber- 
ration free and has NA=0.68, a=0.3 and a reduction factor R=4. Initially, isolated phase 
defects are considered and their ability to transmit the electromagnetic field while shifting 
its phase is determined as a function of their lateral size. Then, the aerial image linewidth 
variation of isolated and dense features with CD^s^lSOnm (IX) is determined in the 
presence of various phase defects. 

In the second part of the Chapter, a domain decomposition strategy appropriate for 
handling wavelength sized phase defects is developed. 

7.1. Phase and brightness of isolated phase defects 

Light that travels the same distance h in two different lossless materials with refrac- 
tive indices n 1? n 2 will have a phase difference in the electric field given by: 

I = ±jz.fi-(n x - n 2 , where X is the illumination wavelength. The sign of this phase differ- 

A 

ence is chosen to accommodate a time retarded field emanating from the high refractive 
index material, where the velocity of propagation is smaller. Therefore, at A.=248nm where 
the refractive index of quartz is n g =1.508 a phase defect (line/post or scratch/void) with 
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height h=240nm should transmit the electric field without changing its amplitude but with 
a phase difference of approximately -180° or +180° with respect to its surrounding. In this 
most frequently assumed ideal geometrical model for phase defects the phase shift depends 
only on the height of the defects and not on their lateral size. In the following, the ideal 1/ 
180° model for phase defects is shown to be incorrect in most cases and the actual trans- 
mission and phase characteristics are seen to depend strongly on the lateral dimensions of 
the defects as well. 

Figure 7-2(a) depicts the normalized amplitude of the near electric field in the vicin- 
ity of an isolated quartz line defect of height h=240nm and width w=200nm, when it is illu- 
minated with a normally incident, linearly polarized plane wave with A=248nm. Figure 7- 
2(b) and (c) show the amplitude and phase of the complex electric field across the observa- 
tion plane. The observation plane is located lOOnm below the quartz line. It is evident that 
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Figure 7-2. Isolated line defect 



(a) Normalized amplitude of the near electric field in the vicinity of an isolated quartz 
line defect of height h=240nm and width w=200nm, when it is illuminated with a nor- 
mally incident, linearly polarized plane wave with A=248nm. (b) and (c) Amplitude 
and phase of the complex electric field across the observation line respectively. 
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the average amplitude across the geometrical shadow of the defect, which is denoted by the 
two vertical lines, is larger than one, immediately drops sharply, and further has shorter 
peaks and valleys. Also, the average phase at the geometrical shadow of this defect is 230°, 
or a change of only -130° in the interval (-180°,180°). Figure 7-3(a) depicts the normalized 
amplitude of the near electric field throughout the cross-section of an isolated quartz 
scratch defect of depth h=240nm and width w=200nm. Figure 7-3(b) and (c) show the 
amplitude and phase of the complex electric field across the observation plane, respec- 
tively. The average light intensity below the defect is much less than the unit intensity of 
the surrounding clear area and the average phase only goes down to 60° instead of 0° as 
expected by geometric optics. There is an even more pronounced deviation from the ideal 
1/180° model in 3D as shown in Figure 7-4 for an isolated quartz post and void of equal 
height (or depth) h=240nm, having a 200nmx200nm square lateral cross-section. The aver- 
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Figure 7-3. Isolated space defect 

(a) Normalized amplitude of the near electric field in the vicinity of an isolated quartz 
scratch defect of depth h=240nm and width w=200nm, when it is illuminated with a 
normally incident, linearly polarized plane wave with X=248nm. (b) and (c) Ampli- 
tude and phase of the complex electric field across the observation line respectively. 
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age brightness (squared amplitude) and phase for the post are 2.1 and -90°, whereas for the 

(b) 



void they are 0.1 and 145° respectively. 

(a) 
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Figure 7-4. Isolated post defect 

(a) Normalized amplitude of the near electric field at the observation plane below a 
200nmx200nm isolated quartz post defect of height h=240nm. (b) Normalized ampli- 
tude of the near electric field at the observation plane below a 200nmx200nm isolated 
void defect of depth h=240nm. (c) and (d) Phase of the complex electric field across 
the cut-lines shown in (a) and (b) respectively. 

The dependence of the average brightness and average phase of phase defects on 
their lateral size is characterized in Figure 7-5 for 2D phase defects (quartz line, quartz 
scratch) and 3D phase defects (quartz post, quartz void). Note that the average phase for 3D 
defects deviates faster from the background phase (0° for posts, lines and 180° for voids, 
scratches) than does the phase for 2D defects. Even for large 2D defects the phase is signif- 
icantly different than the background phase. Figure 7-5(c) shows that the average bright- 
ness of protruding phase defects (posts, lines) is larger than 1 and posts (3D) appear much 
brighter than lines (2D). The average brightness of scratches and voids is seen in Figure 7- 
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5(d) to be less than 1 and it s attains a minimum at a defect size of approximately one free- 
space wavelength. The fact that voids appear darker than posts (or equivalently scratches 
are darker than quartz lines) can be understood by viewing the phase defect as four phase 
edges in close proximity, where light depletes the air-filled regions and concentrates in the 
quartz-filled regions. 
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Figure 7-5. The effect of defect size 

(a) and (b) Average phase and (c) and (d) average brightness (squared amplitude) of 
the complex electric near field underneath a phase defect versus its size for lines/posts 
and scratches/voids. The height (depth) is held constant in all cases at h=240nm. 



The minimum intensity of the aerial image for isolated posts and voids is plotted in 
Figure 7-6 with respect to defect size for phase defects with square lateral shape. Posts 
larger than 150nmxl50nm cause a dip in intensity almost twice as large as the dip caused 

by voids of equal size. The ideal 1/180° model for phase defects, which does not distin- 
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guish between posts and voids, predicts an intensity dip that is between that of a post and a 
void. The minimum intensity of the aerial image for isolated square Chromium defects is 
also plotted in Figure 7-6 versus their size for comparison. The fact that voids smaller than 
150nmxl50nm approach this curve indicates that they are acting almost like opaque defects 
of equal size. 
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Figure 7-6. Minimum intensity of the aerial image for isolated posts and voids 
with respect to defect size for phase defects with square lateral shape and 
h=240nm. The minimum intensity caused by Chromium (opaque) defects is 

also shown for comparison. 



7.2. Simulation results for defective alt. PSMs 

The aerial image degradation of the critical dimension (CD) of 150nm (IX) isolated 
lines and dense line/space patterns was evaluated in the presence of 3D phase defects. 
Quartz posts and voids of various shapes, sizes and heights were placed at various distances 
from the features, while varying the amount of underetch and polarization direction of the 
illumination. The CD aerial was evaluated at the 30% level of the normalized intensity and a 
comparison was made in each case between the worse CD aerial of the defective mask and 
that of the defect-free reference mask. 
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Figure 7-7(a) depicts the fractional change of CD aerial ((CD defective -CD reference )/ 

CD reference) versus defect size. The triangles are the 1/180° geometrical model and serve 
as a point of reference. A 50nmx50nm (IX) post or void is predicted to produce a 10% lin- 
ewidth variation. Void defects are less than 50% as harmful to the aerial image than posts 
of equal size. For example, a 50nmx50nm (IX) void results in a 5% aerial image linewidth 
variation in a dense line/space pattern, whereas a post of the same size results in a poten- 
tially catastrophic 13% ! The ideal 1/180° model for phase defects (triangles) overestimates 
the linewidth variation caused by a void defect by as much as 100% and underestimates the 
linewidth variation caused by a post defect by as much as 40%. These results are in excel- 
lent agreement with the recently presented experimental evidence by Chen et aL [23]. Fully 
buried defects cause minimal linewidth variations and should not be a problem, but when 
they become so large that one of their phase edges exceeds the Chromium edge of the fea- 
ture they can result in significant (>10%) linewidth variations. Additional simulations also 
suggested that the above observations hold regardless of the polarization direction of the 
incident radiation and the exact shape of the defect. In particular, various small shape per- 
turbations were tested and they resulted in almost equal linewidth variations, which renders 
the simplification of parallelepipedly shaped defects valid. This result was independently 
verified by Pistor [67]. 

The fraction of underetch (with respect to the full feature size) was found to have 
no impact on the linewidth variation of a defective mask, unless no underetch is applied, 
in which case the interaction of the phase edge and the defect result in an abrupt increase 
in the linewidth variation. This is shown in Figure 7-7(b), where for a 70nmx70nm defect 
the fractional change of CD aerial jumps from 30% for 50% underetch to 100% (bridging of 
the lines) for no underetch. 

The effect of the defect height on the printability of defects is explored in Figure 7- 
7(c). The fractional change of CD aerial is plotted versus the deviation of the defect phase 

from 1 80° (a defect shorter than h=240nm corresponds to a height-phase smaller than 1 80°, 
therefore the shorter it is the larger its deviation from 180°). The fractional change of 
CDaerial declines slowly when the phase deviation increases, or equivalently when the 
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defect becomes shorter. Similar observations can be found in the recently presented exper- 
iments in [23]. There is an immediate impact of this observation in the repair process, in 
that the partial removal of the phase defect is not as effective as predicted by the 1/180° 



model. 
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Figure 7-7. Defect printability dependence 

(a) Fractional change of CD aerial at the 30% intensity level with respect to defect size 
for 150nm (IX) isolated lines and dense line/space patterns. The impact of posts, 
voids and buried posts under different polarizations (TE or TM) is depicted. The 
effect of ideal 1/180° phase defects (denoted by the triangles) is also shown for com- 
parison. Defect non-planarity and sidewall angle changes produced only small addi- 
tional CD variations as indicated by the clusters of points so labeled, (b)-(d) 
Fractional change of CD aeriaI versus the amount of underetch, phase deviation (defect 
height) and proximity to the feature respectively. 



The fractional change of CD aerial is plotted with respect to the location of the defect 
in Figure 7-7(d). For the dense line/space pattern the location of the defect is irrelevant for 
the linewidth change, whereas for an isolated line placing the defect closer to the line 
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appears to result in a larger linewidth error. This was verified experimentally by Tejnil et 
al in [93]. 

7.3. Efficient phase defect modeling via DDM 

The bright future for alt. PSMs as a major resolution enhancement technique (RET) 
has been plagued primarily by the difficulty to reliably inspect them, flag locations that 
phase defects are present and subsequently repair the defective locations. State-of-the-art 
inspection systems utilize a focused laser beam that can be either at-wavelength (up to 
^=248nm), meaning that the wavelength of the inspection beam is the same as the wave- 
length that the alt. PSM was designed for 1 , or emit light of a larger wavelength. This beam 
scans the entire mask and the reflected or transmitted scattered fields from the particular 
mask location under inspection are imaged with a very high resolution imaging system that 
typically has a reduction factor of R=l and NA close to 1. But here is the caveat: After the 
signal from the inspection tool for a particular mask location is available it needs to be 
compared to something in order to determine whether or not a phase-defect (or other 
defect) is present within that location. Such a comparison is in general very difficult, 
because not only does it need to be accurate in order to reliably flag defective locations, 
but it also needs to be rapid enough for the entire reticle to be inspected in a reasonable 
amount of time. 

The most promising direction that has been adopted by researchers working on the 
problem has been what is known as a die-to-database comparison. By that, it is meant that 
the local signal is compared to the signal that the same mask location would produce if it 
were defect-free. The problem now shifts to building the database of signals from defect- 
free masks. Note that all possible mask geometries that will be encountered in every 
inspected reticle need to exist in the database. One way to building the database would be 
to produce a (set of) test reticle(s) that includes every single geometry situation that is 
anticipated to exist in all designed layouts to come. This (set of) test reticle(s) is then metic- 

I . Recall that the depths of the phase-wells in an alt. PSM depend on the wavelength and on the refractive 
index of the glass substrate, which is in turn a strong function of wavelength. Therefore, although the design 
leads to the required 180° path length difference between alternating apertures "at-wavelength", for differ- 
ent wavelengths this is not the case. 



155 



c 



c 



ulously inspected with a reliable tool (maybe an atomic force microscope - AFM) to guar- 
antee that it is defect-free and finally all signals from the different locations are collected 
and placed in the vault (database). Clearly, the immensity of this task has to do not only 
with the a priori precise anticipation of what layout situations to take provision for, not only 
with manufacturing the test reticle(s), characterizing it and building the database, but also 
with the sheer volume of data that such a database would contain. Moreover, a precise diag- 
nosis of the type, size and location of a phase defect requires an even larger database of 
benchmark signals, where every possible combination of defect type, size and location in 
every possible layout configuration needs to exist in the database. Performing such a task 
in a way similar to the one outlined is impractical. 

Simulation can again come to the rescue, at least in principle. If a simulation tool 
that can rapidly and accurately predict the expected signal resulting from a defect-free or 
defective location on the mask is available, then the task of building the database is simpli- 
fied tremendously. Moreover, if such a tool is really fast (and accurate of course) then a 
database is not even needed! The required benchmark signal from the defect-free layout is 
generated in-situ, while the inspection system gathers measurements. The problem with 
this solution is that to date there exists no simulation tool that possesses both the rapidity 
and accuracy properties. The family of simulation tools that rigorously solve Maxwell's 
equations around the reticle (such as TEMPEST) and subsequently use a vectorial formu- 
lation for the image formation to calculate the expected optical signal from the inspection 
system are impractical because of large memory and time requirements. On the other hand 
the accuracy of speedier simulation programs that circumvent the solution of Maxwell's 
equations around the mask is unacceptable. 

Yet another important consideration of an inspection system for alt. PSMs should 
be the following: Phase-defects that are not critically affecting the image that a projection 
printing tool will produce on the wafer surface should be discarded, i.e. not flagged and not 
repaired. But this is also a difficult task, since the only bullet-proof way of achieving this 
would be to expose wafers with the suspect defect present and examine if the printed resist 
images have intolerable artifacts. Clearly, such a process now involves the close coopera- 
tion of the inspection system with the projection printing tool. However, an accurate and 
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properly calibrated simulation tool can alternatively assess the tendency of the defect to 
print or cause other image artifacts. Consequently, the costly step of exposing wafers can 
be avoided and a decision of repairing or not the defect rests fully upon the simulation 
result. Again, to date, results from no simulation tool, that is fast enough in calculating 
expected images with phase defects present, can be trusted to base repair decisions on. 

The domain decomposition methods developed in Chapters 4 - 6 are appropriate for 
solving the problem of the rapid and accurate evaluation of the benchmark signal from non- 
defective masks for the inspection system. In this Section, a domain decomposition method 
is presented that addresses the problem of the defect printability assessment, after a defect 
is previously found. The generation of benchmark signals from defective mask locations 
can also be based on the domain decomposition technique, although more work will be 
needed to adapt it and test it for the high-NA inspection optics with R=l. As will soon be 
obvious, the speed and efficiency of the method are inherited by virtue of the domain 
decomposition methods presented in Chapters 4-6, although for simplicity in the present 
discussion the qr-DDM and edge-DDM of Chapters 4 and 6 are not explicitly employed in 
the following. In the simulation examples of this Section the illumination wavelength is 
193nm and the phase-wells are designed accordingly, so as to provide the required phase 
shifts for that wavelength. 

Consider the example depicted in Figure 7-8. In (a) the near scattered field below 

a small layout of a 0°/180° alt. PSM is shown when a square 120nm x 120nm (4X), 180° 
post phase defect is present at various locations within the layout. Note that although nei- 
ther the layout nor the phase defect geometry change, a different simulation altogether is 
needed because of the different relative positions of the defect within the layout. Each one 
of the nine 3D simulations in the matrix shown in (a) take approximately lOhrs on a 
450MHz CPU and utilize 300Mb of memory. Subsequently, the aerial images for all nine 
situations are calculated for an imaging system with R=4, NA=0.75 and a=0.3. Then, the 
fractional change of CD aeriai at the middle and at the end of the line because of the defect 
is found for all nine defect locations and is shown in Figure 7-8(b). From the plot in 
Figure 7-8(b) it can be deduced that when the defect is located in the middle of the phase- 
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wells or is mostly "tucked" underneath the absorbing layer (columns 3 and 1 of Figure 7- 
8(a) respectively) the CD aerial is not critically affected, but locations such as those in 
column 2 of Figure 7-8(a) cause larger CD aeria , variations. In the former cases the defects 
could probably be discarded whereas in the latter repair seems necessary. Clearly, the 
luxury of the lOhr-long simulations is only available as a proof of concept and such an 
approach is not viable for full-chip characterizations. 

(a) (b) 
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Figure 7-8. Defect printability assessment (the laborious way) 

A square 120nm x 120nm phase defect (post) with 180° height causes different CD 
variations depending on its relative location with respect to the center line. In general, 
larger CD variation is evidenced when the defect is close to the line (column 2 in (a)).' 
However, when the defect is partially covered by the Cr-layer of the line, smaller CD 
changes are observed. This is quantified in (b). The imaging system parameters are: 
X=193nm, NA=0.75, R=4, a=0.3. Note that each simulation in (a) requires ~10hrs on 
a 450MHz CPU. 



Now, consider the paradigm shown in Figure 7-9. A domain decomposition method 
is again invoked, where instead of a one-step simulation of the defective layout two sepa- 
rate rigorous simulations are performed, one with just the defect-free layout and one with 
just the defect in clear surroundings. The scattered field below the defect is then shifted 
accordingly, such that the defect is effectively placed at the location it appears in the defec- 
tive mask, and the uniform background light is subtracted, so as to retrieve the (complex) 
signal of the local perturbation that the defect causes. This perturbation is added to the 
defect-free layout signal for an approximation to the signal from the rigorous simulation of 
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the defective mask. The non-defective mask simulation can be performed rapidly using the 
decomposition methods of Chapters 4-6. 




Figure 7-9. Domain decomposition method when a defect is present 

A strategy can now be devised, where a database of the near electromagnetic field 
scattered from isolated phase-defects in clear surroundings is created for all possible poten- 
tial defects that are critical. The volume of such a database can be reduced based on the 
earlier observations that defects with similar footprints (lateral dimensions) and small 
shape perturbations are equivalent from an electromagnetic point-of-view. 

The suggested decomposition method will be introduced through the following 
example: Suppose a 150nm (IX), 1:1.5 semi-dense, 9G°/270° alt. PSM contact mask. The 
layout and the scattered field across the observation plane under normal incidence E y (TE) 
polarization are shown in Figure 7-10. 
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Figure 7-10. Example: 150nm (IX), 1:1.5 semi-dense, 90°/270° alt. PSM contact mask 
(a) Layout and (b) amplitude of scattered field below the mask at the observation plane. 
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Now, suppose that a 200nm by 200nm (4X) 90° post defect is present in the center 
of the bottom left hole. According to the paradigm of Figure 7-9, the mask scattering sim- 
ulation can be broken up into two constituent parts. This approach is shown in Figure 7-11. 
In (a) the phase defect is simulated in an isolated configuration at the location that it appears 
in the layout. In (b) and (c) the near field across the observation plane below the defective 
alt. PSM is shown, using the domain decomposition method and the complete mask respec- 
tively, The differences are not discernible from the plots of (b) and (c) and one has to look 
at the amplitude of the error in (d). Again, not only is the error level of the decomposition 
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Figure 7-11. A 90° phase defect present in the center of the bottom-left hole 

Amplitude of scattered field across the observation plane for (a) isolated phase defect 
(b) defective mask simulated with the DDM and (c) rigorous (one-step) simulation of 
defective mask. The amplitude of the error of the DDM compared to the rigorous sim- 
ulation ,s shown in (d). The normalized mean square error calculated at the bottom 
left (defective) quadrant is only 0.7%. 
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low, but the frequency content is very high, so that it will be filtered out by the imaging 
system. 

Next, suppose that the same 90° post defect is present in the center of the bottom 
right hole. The important observation here is that no new simulation is necessary for the 
decomposition method* The scattering of the isolated phase defect can be recycled from the 
database, i.e. from Figure 7-1 1(a). However, the scattered field is now properly shifted in 
the lateral x-direction and also (a key step) it is propagated in the z-direction according to 
the discussion in Appendix A. This results in the plot of Figure 7-12(a), where a larger dif- 

(a) (b) 





(d) 




Figure 7-12. A 90° phase defect present in the center of the bottom-right hole 

Amplitude of scattered field across the observation plane for (a) isolated phase 
defect, (b) defective mask simulated with the DDM and (c) rigorous (one-step) 
simulation of defective mask. Note that the scattered field in (a) does not result 
from a new simulation. It comes from the field plot calculated in Figure 7-1 1(a) 
that is properly shifted in x and propagated in z to account for the different loca- 
tion of the defect. The amplitude of the error of the DDM compared to the rigor- 
ous simulation is shown in (d). The normalized mean square error calculated at 
the bottom right (defective) quadrant is only 0.5%. 
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fraction spreading is evident. In (b) and (c) the near field across the observation plane below 
the defective alt. PSM is shown, using the domain decomposition method and the complete 
mask respectively. The amplitude of the error is shown in (d) and similar observations 
apply in Figure 7-1 1(d). 

The near fields across cut-lines passing through the center of the bottom two holes 
are shown for the two defective masks in Figure 7-13, comparing the rigorous (r-mask) and 
the decomposition methods (DDM). The non-defective case is also plotted for comparison. 
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Figure 7-13. Comparison of the decomposition method with the rigorous solution 

(a) Phase defect in the 90° hole and (b) phase defect in the 270° hole. The cut-lines 
are across x, passing through the center of the defective holes. 



Next, the images of both rigorous and DDM approaches, for both defective masks 
are compared in Figure 7-14. The optical system parameters are A=193nm, R=4, NA=0.75 
and a=0.3. The agreement between the r-mask and the DDM is good, a normalized mean 
square error of less than 0.3% is incurred (measured only in the area of the defective hole). 
The k-mask model is insufficient for capturing these effects and is only shown for compar- 



ison. 
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Figure 7-14. Aerial images of defective masks of Figure 7-11 and Figure 7-12 

(a) Phase defect in the 90° hole (bottom-left) and (b) phase defect in the 270° hole 
(bottom-right). Imaging system parameters X=193nm, R=4, NA=0.75, a=0.3. 



The reason of the excellent agreement in the predictions of the aerial image can be 
traced back to the plots of Figure 7- 11(d) and Figure 7- 12(d). The high spatial frequency 
variation of the error means that it is mapped at the extremities of the spectrum of propa- 
gating plane waves and is not collected by the optical system. The normalized spectra 
(magnitude) of the near field errors are shown in Figure 7-15. Observe that the implied 
accuracy from these plots should be enough even for R=l and NA-0.8-0.9, which would 
be sufficient for inspection simulations. 
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Figure 7-15. Spectra plots (magnitude) of the errors of Figure 7-ll(d) and Figure 7-12(d) 

The error is concentrated at the extremities of the spectrum of propagating plane waves 
indicated by the outer circle and does not exceed 6%. However, the projection printing 
system (with A=193nm, R=4, NA=0.75, a=0.3) uses only the diffraction orders indicated 
by the inner circle, where the error levels are minimal. 

7.4. Conclusions 

Rigorous 3D electromagnetic simulation is suitable for understanding phase defects 
and their impact on aerial images of features. Isolated phase defects were found to deviate 
significantly from the ideal geometric-optics model in which the fields are transmitted with 
unaltered amplitude and 180° shift in phase. The dependence of the amplitude and phase of 
the transmitted field on the lateral size of defects for both isolated quartz posts and voids 
was found. Small voids behave more like opaque defects of equal size, whereas posts can 
result in a twice as large intensity drop as voids. This results in posts being more than twice 
as harmful to the aerial image than voids of equal size, when phase defects are located close 
to features. In particular, a 50nmx50nm (IX) void results in an acceptable 5% linewidth 
change whereas a post of the same size causes an unacceptable 13% variation. Fully buried 
defects result in minor linewidth variations. Light polarization conditions, the amount of 
underetch and the defect height are 2nd order important in defect printability. The linewidth 
variations of the aerial images were also insensitive to the exact defect shape, which vali- 
dates the use of parallelepipeds as a sufficient shape to address defect printability issues. 
The interesting part of this research work is that most of these results were first observed 
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through simulations and recently experimental evidence is becoming available that vali- 
dates the above observations. 

In the second part of the Chapter a suitable extension of the domain decomposition 
framework for the simulation of alternating phase shift mask with phase defects was devel- 
oped. The value of this approach is believed to be significant in the problem of rapid 
assessment of defect printability once a defect has been found, or in the die-to-database 
comparisons of inspection systems when deciding about the existence or not of a phase- 
defect. 
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AS integrated circuit technology marches towards a minimum feature resolution of 
less than lOOnm and k, factors well below 0.5, the already important role of process simu- 
lation and modeling of the various physical phenomena that take place is expected to grow. 
In particular, many of the proposed resolution enhancement techniques (RET), such as 
optical proximity correction (OPC), phase shifting masks (PSM), concurrent illumination/ 
mask optimization, etc., that are expected to prolong the life of optical photolithography, 
rely heavily on modeling and computer simulations. On one hand, rigorous methods are 
based on the solution of Maxwell's equations for the electromagnetic field in three-dimen- 
sional structures and possess excellent accuracy but their computational requirements 
(memory and speed limitations) render them inappropriate for large scale problems. On the 
other hand, the simplifying assumptions made for the mask scattering in scalar or vector 
diffraction theory formulations are inaccurate for a number of problems of practical inter- 
est. A central part of this work has been an attempt to bridge the gap between rigorous sim- 
ulation methods and simpler methods based on diffraction theory for the photomask 
imaging process. 

In this final Chapter, a summary of the most important results of this research is 
given first, and directions for future research and further development of the presented sim- 
ulation methods are discussed next. 

8.1. Summary of important results 

The most important contribution of this work is doubtlessly the domain decompo- 
sition methods (DDMs) developed and used in Chapters 4 through 6 and partly in Chapter 
7. Two central ideas were sequentially exploited for the development of DDMs: First, orig- 
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inating from the linearity of the Kirchhoff-Fresnel diffraction integral, a division of a mul- 
tiple-opening diffracting screen (photomask with a number of features carved on it) was 
performed and the problem was broken up in elemental, single-opening pieces. Each one 
of these pieces was treated separately and it was shown that the near scattered field below 
such a single-opening mask, that would normally require solving Maxwell's equations on 
a three-dimensional domain, can be accurately approximated by combining the results of 
the scattered field from two two-dimensional structures. The errors incurred near corners 
and the fact that light polarizations other than the excited one were neglected were shown 
to be insignificant. The reason for that was traced to the energy concentration of both 
corner effects and perpendicular polarizations into the higher spatial frequencies of the dif- 
fraction spectra. The role of the imaging system (with NA<1 and R~4 or 5) as a filter for 
the higher spatial frequencies is such that these errors are filtered-out and they do not con- 
found the simulation of image formation. The second idea had to do with the use of simple 
functions (like the rect-function, raised-cosine, or Gaussian) that substitute the true, rigor- 
ously calculated, diffraction field profiles below the elemental two-dimensional geome- 
tries in a way that they perfectly match (to within some level of accuracy that was shown 
to be better than 99%) their through-the-lens spectra. A good spectral match within the part 
of the spectrum that is actually used by the projection system for the image formation is all 
that is required for these models to work. And so they do. 

Being able to handle arbitrary two-dimensional layouts is key for the practical 
aspect of domain decomposition methods. To this extent a methodology was developed in 
Chapter 4 that involved breaking up arbitrary (Manhattan-type) polygons into the set of 
maximally overlapping rectangles and it was explained why this break-up method is supe- 
rior from an electromagnetic point-of-view to the one where the polygon is divided into a 
set of mutually disjoint rectangles. 

Cross-talk effects between neighboring apertures in alternating phase shift masks 
were the theme of Chapter 5. A successful model for their prediction was build based upon 

a single simulation for the electromagnetic field scattered from a 90° corner at a glass/air 
interface, from which it became clear that the edges at the bottom of the phase-wells are 
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the primary sources of cross-talk. Two paths, of almost equal strength, that cause cross-talk 
were identified and the amount of cross-talk was quantified for various mask geometries, 
edge profiles and field polarizations. Recipes for cross-talk elimination were also sug- 
gested. Proper modeling of cross-talk effects was shown to be an integral part of domain 
decomposition methods. In particular, examples of two-dimensional mask layouts where 
cross-talk cannot be neglected were presented. 

In Chapter 6 the edge domain decomposition method was developed. The main 
ideas were similar to those in Chapter 4, but here the decomposition was based on edges, 
i.e. the electromagnetic field diffraction profiles from isolated edges. This edge-decompo- 
sition allowed for tremendous versatility, since only a small number of two-dimensional 
simulations of edges is enough for the reconstruction of accurate approximations of the 
near field below two-dimensional arbitrary mask layouts. A 12|im by 16|Lim layout of a 
three level alternating PSM that is out of reach for fully rigorous methods was simulated in 
under lmin. Matching the spectra of the fields from edge-diffractions with simple piece- 
wise constant models was also possible here and further accelerated the calculations. 
Although in essence the method is not different than the domain decomposition method of 
Chapter 4, its algorithmic implementation for arbitrary 2D layouts is easier and it requires 
far less rigorous simulations to build the library of required edges, as opposed to building 
a library for every possible feature size required by the decomposition method in Chapter 4. 

Useful design strategies for optical proximity correction can be found in Chapter 3. 
A critical comparison between rigorous simulation methods and scalar methods based on 
the Kirchhoff approximation for the mask was carried out and it was found that the -10% 
bias in the results for binary Cr-on-glass (COG) masks obtained from the two different 
methods can be adequately explained by the use of a perturbation model that incorporates 
the electromagnetic behavior of light penetrating through small (wavelength-sized) holes. 

Differences in the two types of defects (posts and voids) encountered in alternating 
phase shift masks were exposed in Chapter 7. The surprising behavior of light passing 
through phase-defects, where the change in phase depends strongly of the footprint (lateral 
size) of the defect for small defects was characterized. This behavior was further used to 
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explain the changes in critical dimension (CD) that phase defects cause depending on their 
geometry and proximity to layout features. Recent experimental results are in excellent 
agreement with those observations. Also, in Chapter 7 a methodology based on domain 
decomposition was proposed for the rapid simulation of defective masks. This methodol- 
ogy is expected to have a significant impact in die-to-database comparisons in mask 
inspection systems and in the rapid evaluation of defect printability. 

8.2. Future research 

Speed-up factors in excess of 200, without compromising the accuracy, were pos- 
sible with the domain decomposition methods of Chapter 4. However, it was stressed that 
a possible algorithmic implementation could involve the pre-calculation of all different 
feature sizes that will be encountered on a layout and creation of a database with their dif- 
fraction characteristics (parameters of the compact mk-mask model). Such an approach is 
expected to lead to a much larger acceleration. The speed-up factor of 172,800 (lsec. 
instead of 2days) achieved in the simulation of a simple (isolated hole and isolated island) 
3|am by 3\xm layout in the case of the edge-DDM should be possible for the DDM in Chap- 
ter 4 also. However, the automation of this task requires considerable effort. Algorithms 
for the rapid synthesis of the diffraction characteristics of arbitrary, complex layouts will 
be needed, especially in order to decide when cross-talk effects are important and automat- 
ically include them. 

Cross-talk effects were not treated in the case of the edge-DDM of Chapter 6. It 
should be obvious that although this method has an easier algorithmic implementation 
compared to the decomposition methods of Chapter 4, including cross-talk effects raises 
the complications to a similar level as in the methods of Chapters 4 and 5. It was also found 
in Chapter 6 that the accuracy of the edge-DDM drops sharply for sub-wavelength sized 
features. Clearly, in that regime a new cross-talk effect emerges, that between the closely 
packed edges. It would be worth to explore if including this cross-talk effect is possible, 
and if so, determine whether the revised edge-DDM is more or less efficient in large scale 
computations than simulating every feature size separately without decomposing it in 
edges, as in the methods of Chapters 4 and 5. 
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There is nothing that is restrictive in the use of the domain decomposition methods 
for normal incidence illumination. This fact was validated with off-axis illumination exam- 
ples, where the methods still perform excellent. However, a considerable effort is required 
in order for the method to become available for off-axis partially coherent illumination, 
when the Hopkins' approximation of constant diffraction orders with respect to the angle 
of incidence breaks down. Such an extension of Hopkins' approximation will be needed for 
example, if the DDM is used in the simulation of images from inspection systems (R=l, 
NA~1). 
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A Angular spectrum decomposition of 
^^electromagnetic waves 



If the complex field distribution of a monochromatic disturbance is Fourier-ana- 
lyzed across any plane, the various spatial Fourier components can be identified as plane 
waves traveling in different directions away from that plane. The field amplitude at any 
other point or across any other parallel plane can be calculated by adding the contributions 
of these plane waves, taking into account the phase shifts they have undergone during 
propagation. The following discussion is very helpful in understanding the physical mean- 
ing of the spectrum (Fourier transform) of the scattered complex electromagnetic field 
across the observation plane below photomasks, and also how the choice of the location of 
the observation plane affects the spectrum, as seen throughout this work. 

A.l. The angular spectrum and its physical interpretation 

Suppose the presence of an electromagnetic wave traveling with a component of 
propagation in the positive z-direction and observed on a transverse (x,y) plane. Let 
U(x,y,0) be the complex field across the z=0 plane. U can represent either one of the three 
cartesian components of the electric or magnetic field. An expression is sought that relates 
the field U(x,y,z) appearing across a plane z>0, parallel to z=0 plane, with the field 
U(x,y,0) across the z=0 plane. 

Across the z=0 plane, the function U has a two-dimensional Fourier transform 
given by: 

eo 

Equation A-1. A(k x , k y ;0) = fj U(x, y, 0)exp[-j(k x x + k y y)]dxdy 
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where k x -27tf x and k y =27tf y . U(x,y,0) can also be written as an inverse Fourier transform 
of its spectrum A(k x ,k y ;0): 



Equation A-2. U(x 9 y, 0) = J J A(*„ ^;0)exp[/(^ + k y)]dk x dk 



To give physical meaning to the functions 
in the integrand of the above integral, consider the 
form of a simple plane wave propagating with 

wave vector k , where k has magnitude 2n/X and 
has direction cosines (a,p,y), as illustrated in 
Figure A-l. Such a plane wave has a complex rep- 
resentation of the form: 




k-vector 



cos^oc y 



Figure A-l. Wave vector k and 
direction cosines (a,P,y) 



Equation A-3. 



> > 



™(x,y,z;t) = exp [/(*• r- tor)] 
where r = xx + yy + zz is a position vector, co is the angular frequency of the monochro- 



matic wave and the wavevector k is: 



2n, 



Equation A-4. k = kjc + fc y y + * z z = y (ax + 0y + y z ) 

Dropping the time dependence, the complex phasor amplitude of the plane wave 
across a constant z-plane is: 

Equation A-5. P(x, y, z) = exp(jk ■ r) = exp[}y (cu + py)J • exp(/y yz) 

Note that the direction cosines are not independent but they are related through: 



Equation A-6. y = Ji . a 2 - p 2 

Therefore across the plane z=0, a complex-exponential function exp + ^y)] 
may be viewed as a unit-amplitude plane wave propagating with direction cosines: 
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In the Fourier decomposition of U, the complex amplitude of the plane-wave com- 
ponent with spatial frequencies (k x ,k y ) is simply A(k x ,k y ;0)dk x dk y , evaluated at (k x =27tcc/ 
X, k y =27iP/A). For this reason, the function A(k x ,k y ;0) (spectrum of U(x,y,0)) expressed 
explicitly as function of a, P is called the angular spectrum of the disturbance U(x,y,0): 

oo 

Equation A-8. A^jj-, $ ;0j = J J U(x, y, 0)exp ["7'2jc^jc + |y)]<My 

-oo 

A.2. Propagation of the spectrum 

Let A(k x ,k y ;z) be the spectrum of the disturbance U(x,y,z) across a plane parallel 
to the (x,y) plane but at a distance z from it: 

oo 

Equation A-9. A{k x , k y ;z) = J J U(x, y, z)e\p[-j(k x x + k y y)]dxdy 

-oo 

If the relation between A(k x ,k y ;0) and A(k x ,k y ;z) can be found, then the effects of 
wave propagation on the spectrum of the disturbance will be evident. Expressing U(x,y,z) 
as the inverse Fourier-transform of its spectrum A(k x ,k y ;z): 

oo 

Equation A-10. U(x, y,z) = j J A(k x , k^cxpljik^ + k y y)]dk x dk y 

-oo 

and considering the fact that U must satisfy the Helmholtz equation at all source-free 
points: 

Equation A-l 1 . V 2 U + k 2 JJ = 0 (Helmholtz eq.) 

it follows that A must satisfy the following differential equation: 

d 2 2 
Equation A-12. —A(k x , k y ;z) + (k Q - k x - k y )A(k x , k y ;z) = 0 

dz 
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where Y^Znfk is the magnitude of the wavevector in the ambient. An elementary solution 
of Equation A-l 2 can be written as: 



Equation A-13. A(k x , k y ;z) = A(k x , k y ,0)cxp(jJk 2 Q - k\ - k] - z ) 

Based on this result, when k\ + k) < k\ = (^j* or equivalently when the direc- 
tion cosines (a,0) satisfy a 2 + 0 2 < 1 , as all true direction cosines must, the effect of prop- 
agation over distance z is simply a change of the relative phases of the various components 
of the angular spectrum. Since each plane-wave component propagates at a different angle, 

w *s u ui^iuutc uciwccn iwo parallel planes, and relative phase delays are 

therefore introduced. 

2 2 2 

However, when k x + k y > k Q or equivalently when the direction cosines (cc,0) sat- 

2 2 

isfy a + p > 1 , a different interpretation is required, a and 0 can no longer be viewed as 
direction cosines. The square root in Equation A-13 is imaginary, and that equation can be 



rewritten as: 



Equation A-14. A {k x , k y ,z) = A{k x , ^;0)exp(- Hz ) 

I 2 2 2 

where M- = *Jk x + k y - k 0 . Since m is a positive real number, these wave components are 
rapidly attenuated by the propagation phenomenon. Such components are called evanes- 
cent waves and they carry no energy away from the transverse plane. 

A.3. Choosing the location of the observation plane 

Throughout this work whenever the scattered field below a photomask is rigorously 
calculated (with TEMPEST) the choice of the observation plane has seemed somewhat 
arbitrary. Although the observation plane has been typically chosen at a distance d = A/4 
(for example at d=50nm for A=193nm in Chapters 4-6) it has been suppressed how the 
choice of the exact location of the observation plane affects the subsequent results. The fol- 
lowing example demonstrates how the spectrum (diffraction orders for periodic masks) of 
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the scattered field changes as the observation plane is moved at different distances from 
the mask, according to the theory presented in Section A.l and Section A.2. 



An isolated 400nm (4X) mask opening with 180° phase shift and 50nm of isotropic 
underetch is shown in Figure A-2(a). The amplitude of the complex electromagnetic field 
being established when a plane wave of A,=193nm is normally incident on the mask is also 
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Figure A-2. Location of observation plane and diffraction spectrum 

An isolated 400nm (4X) mask opening with 180° phase shift and 50nm of isotropic 
underetch is shown in (a) with the amplitude of the field everywhere in the simulation 
domain under a normally incident plane wave excitation at A=193nm. Although the 
profile of the scattered near field depends strongly on the exact choice of the location 
d of the observation plane, as shown in (b), the amplitude of the spectrum does not 
depend on d, as shown in (c). However, the phase of the higher diffraction orders is 
seen to vary rapidly with d, whereas the phase of the lower diffraction orders is quite 
insensitive to d. These results are in exact agreement with Equation A- 13. 
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shown in Figure A-2(a) everywhere in the simulation domain. Four observation planes sep- 
arated by 50nm (-A/4) are also indicated, where the one nearest to the mask is located 25nm 
below the Cr-based absorption layer. There are visible changes in the amplitude of the scat- 
tered field at the different observation planes, as shown in Figure A-2(b). This however is 
not significant as one might expect from the previous treatment on the propagation of the 
spectrum. Consistent with Equation A- 13, the amplitude of the spectrum does not change 
as the distance d of the observation plane from the mask changes. This is depicted in 
Figure A-2(c). The phase of the diffraction orders shown in Figure A-2(d) is seen to depend 
on the location of the observation plane. Again, this is consistent with Equation A- 13 and 
the phase of the higher orders is affected more from the shift in the location of the observa- 
tion plane than the phase of the lower orders. Owing to this observation and the fact that a 
typical projection imaging system for semiconductor manufacturing has a reduction factor 
of at least 4X it turns out that the image will be quite insensitive on the exact location of the 
observation plane (provided that this plane is kept nearer than approximately one wave- 
length from the mask). 
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